degrees <-> radians conversion clean-up

src/chart/deltaE.c

@@ -1,6 +1,6 @@
 /*
  *    This file is part of darktable,
- *    Copyright (C) 2016-2020 darktable developers.
+ *    Copyright (C) 2016-2025 darktable developers.
  *
  *    darktable is free software: you can redistribute it and/or modify
  *    it under the terms of the GNU General Public License as published by
@@ -17,24 +17,21 @@
  */
 
 #include "deltaE.h"
+#include "common/math.h"
 
 #ifndef _XOPEN_SOURCE
 #define _XOPEN_SOURCE 600
 #endif
 
 #include <math.h>
 
-
-#define DEG2RAD(deg) (deg * M_PI / 180.0)
-#define RAD2DEG(rad) (rad * 180.0 / M_PI)
-
 // http://www.brucelindbloom.com/index.html?Eqn_DeltaE_CIE76.html
 float dt_colorspaces_deltaE_1976(dt_aligned_pixel_t Lab0, dt_aligned_pixel_t Lab1)
 {
   float dE = 0.0;
   for(int i = 0; i < 3; i++)
   {
-    float difference = Lab0[i] - Lab1[i];
+    const float difference = Lab0[i] - Lab1[i];
     dE += difference * difference;
   }
   return sqrtf(dE);
@@ -43,27 +40,27 @@ float dt_colorspaces_deltaE_1976(dt_aligned_pixel_t Lab0, dt_aligned_pixel_t Lab
 // http://www.brucelindbloom.com/index.html?Eqn_DeltaE_CIE2000.html
 float dt_colorspaces_deltaE_2000(dt_aligned_pixel_t Lab0, dt_aligned_pixel_t Lab1)
 {
-  float L_ip = (Lab0[0] + Lab1[0]) * 0.5;
-  float C1 = sqrtf(Lab0[1] * Lab0[1] + Lab0[2] * Lab0[2]);
-  float C2 = sqrtf(Lab1[1] * Lab1[1] + Lab1[2] * Lab1[2]);
-  float C_i = (C1 + C2) * 0.5;
-  float G = (1.0 - sqrtf(powf(C_i, 7) / (powf(C_i, 7) + powf(25, 7)))) * 0.5;
-  float a1_p = Lab0[1] * (1 + G);
-  float a2_p = Lab1[1] * (1 + G);
-  float C1_p = sqrtf(a1_p * a1_p + Lab0[2] * Lab0[2]);
-  float C2_p = sqrtf(a2_p * a2_p + Lab1[2] * Lab1[2]);
-  float C_ip = (C1_p + C2_p) * 0.5;
-  float h1_p = RAD2DEG(atan2f(Lab0[2], a1_p));
+  const float L_ip = (Lab0[0] + Lab1[0]) * 0.5;
+  const float C1 = sqrtf(Lab0[1] * Lab0[1] + Lab0[2] * Lab0[2]);
+  const float C2 = sqrtf(Lab1[1] * Lab1[1] + Lab1[2] * Lab1[2]);
+  const float C_i = (C1 + C2) * 0.5;
+  const float G = (1.0 - sqrtf(powf(C_i, 7) / (powf(C_i, 7) + powf(25, 7)))) * 0.5;
+  const float a1_p = Lab0[1] * (1 + G);
+  const float a2_p = Lab1[1] * (1 + G);
+  const float C1_p = sqrtf(a1_p * a1_p + Lab0[2] * Lab0[2]);
+  const float C2_p = sqrtf(a2_p * a2_p + Lab1[2] * Lab1[2]);
+  const float C_ip = (C1_p + C2_p) * 0.5;
+  float h1_p = rad2degf(atan2f(Lab0[2], a1_p));
   if(h1_p < 0) h1_p += 360.0;
-  float h2_p = RAD2DEG(atan2f(Lab1[2], a2_p));
+  float h2_p = rad2degf(atan2f(Lab1[2], a2_p));
   if(h2_p < 0) h2_p += 360.0;
   float H_ip;
   if(fabsf(h1_p - h2_p) > 180.0)
     H_ip = (h1_p + h2_p + 360.0) * 0.5;
   else
     H_ip = (h1_p + h2_p) * 0.5;
-  float T = 1.0 - 0.17 * cosf(DEG2RAD(H_ip - 30.0)) + 0.24 * cosf(DEG2RAD(2.0 * H_ip))
-            + 0.32 * cosf(DEG2RAD(3.0 * H_ip + 6.0)) - 0.20 * cosf(DEG2RAD(4.0 * H_ip - 63.0));
+  const float T = 1.f - 0.17f * cosf(deg2radf(H_ip - 30.f)) + 0.24f * cosf(deg2radf(2.f * H_ip))
+            + 0.32 * cosf(deg2radf(3.f * H_ip + 6.f)) - 0.2f * cosf(deg2radf(4.f * H_ip - 63.f));
   float dh_p = h2_p - h1_p;
   if(fabsf(dh_p) > 180.0)
   {
@@ -72,20 +69,20 @@ float dt_colorspaces_deltaE_2000(dt_aligned_pixel_t Lab0, dt_aligned_pixel_t Lab
     else
       dh_p -= 360.0;
   }
-  float dL_p = Lab1[0] - Lab0[0];
-  float dC_p = C2_p - C1_p;
-  float dH_p = 2.0 * sqrtf(C1_p * C2_p) * sinf(DEG2RAD(dh_p * 0.5));
-  float SL = 1.0 + ((0.015 * (L_ip - 50.0) * (L_ip - 50.0)) / sqrtf(20.0 + (L_ip - 50.0) * (L_ip - 50.0)));
-  float SC = 1.0 + 0.045 * C_ip;
-  float SH = 1.0 + 0.015 * C_ip * T;
-  float dtheta = 30.0 * expf(-1.0 * ((H_ip - 275.0) / 25.0) * ((H_ip - 275.0) / 25.0));
-  float RC = 2.0 * sqrtf(powf(C_ip, 7) / (powf(C_ip, 7) + powf(25, 7)));
-  float RT = -1.0 * RC * sinf(DEG2RAD(2.0 * dtheta));
-  float KL = 1.0;
-  float KC = 1.0;
-  float KH = 1.0;
+  const float dL_p = Lab1[0] - Lab0[0];
+  const float dC_p = C2_p - C1_p;
+  const float dH_p = 2.0 * sqrtf(C1_p * C2_p) * sinf(deg2radf(dh_p * 0.5f));
+  const float SL = 1.0 + ((0.015 * (L_ip - 50.0) * (L_ip - 50.0)) / sqrtf(20.0 + (L_ip - 50.0) * (L_ip - 50.0)));
+  const float SC = 1.0 + 0.045 * C_ip;
+  const float SH = 1.0 + 0.015 * C_ip * T;
+  const float dtheta = 30.0 * expf(-1.0 * ((H_ip - 275.0) / 25.0) * ((H_ip - 275.0) / 25.0));
+  const float RC = 2.0 * sqrtf(powf(C_ip, 7) / (powf(C_ip, 7) + powf(25, 7)));
+  const float RT = -1.0 * RC * sinf(deg2radf(2.f * dtheta));
+  const float KL = 1.0;
+  const float KC = 1.0;
+  const float KH = 1.0;
 
-  float dE = sqrtf((dL_p / (KL * SL)) * (dL_p / (KL * SL)) + (dC_p / (KC * SC)) * (dC_p / (KC * SC))
+  const float dE = sqrtf((dL_p / (KL * SL)) * (dL_p / (KL * SL)) + (dC_p / (KC * SC)) * (dC_p / (KC * SC))
                    + (dH_p / (KH * SH)) * (dH_p / (KH * SH)) + RT * (dC_p / (KC * SC)) * (dH_p / (KH * SH)));
   return dE;
 }

src/chart/thinplate.c

@@ -1,6 +1,6 @@
 /*
  *    This file is part of darktable,
- *    Copyright (C) 2016-2020 darktable developers.
+ *    Copyright (C) 2016-2025 darktable developers.
  *
  *    darktable is free software: you can redistribute it and/or modify
  *    it under the terms of the GNU General Public License as published by
@@ -25,6 +25,7 @@
 
 #include "chart/thinplate.h"
 #include "chart/deltaE.h"
+#include "common/math.h"
 #include "iop/svd.h"
 
 #include <assert.h>
@@ -36,9 +37,9 @@
 // #define EXACT       // use full solve instead of dot in inner loop
 
 // very fast, very approximate
-static inline float __attribute__((__unused__)) fasterlog(float x)
+static inline float __attribute__((__unused__)) fasterlog(const float x)
 {
-  union { float f; uint32_t i; } vx = { x };
+  const union { float f; uint32_t i; } vx = { x };
   float y = vx.i;
   y *= 8.2629582881927490e-8f;
   return y - 87.989971088f;
@@ -110,7 +111,7 @@ static inline double compute_error(
   return err;
 }
 
-static inline int solve(double *As, double *w, double *v, const double *b, double *coeff, int wd, int s, int S)
+static inline int solve(double *As, double *w, double *v, const double *b, double *coeff, const int wd, const int s, const int S)
 {
   // A'[wd][s+1] = u[wd][s+1] diag(w[s+1]) v[s+1][s+1]^t
   //
@@ -432,12 +433,10 @@ int thinplate_match(const tonecurve_t *curve, // tonecurve to apply after this (
 
 #pragma GCC diagnostic pop
 
-float thinplate_color_pos(float L, float a, float b)
+float thinplate_color_pos(const float L, const float a, const float b)
 {
-  const float pi = 3.14153f; // clearly true.
-  const float h = atan2f(b, a) + pi;
-  // const float C = sqrtf(a*a + b*b);
-  const int sector = 4.0f * h / (2.0f * pi);
+  const float h = atan2f(b, a) + M_PI_F;
+  const int sector = 4.0f * h / (2.0f * M_PI_F);
   return 256.0 * sector + L; // C;
 }
 

src/common/gpx.c

@@ -1,6 +1,6 @@
 /*
     This file is part of darktable,
-    Copyright (C) 2011-2023 darktable developers.
+    Copyright (C) 2011-2025 darktable developers.
 
     darktable is free software: you can redistribute it and/or modify
     it under the terms of the GNU General Public License as published by
@@ -61,14 +61,14 @@ static GMarkupParser _gpx_parser
     = { _gpx_parser_start_element, _gpx_parser_end_element, _gpx_parser_text, NULL, NULL };
 
 
-static gint _sort_track(gconstpointer a, gconstpointer b)
+static gint _sort_track(const gconstpointer a, const gconstpointer b)
 {
   const dt_gpx_track_point_t *pa = (const dt_gpx_track_point_t *)a;
   const dt_gpx_track_point_t *pb = (const dt_gpx_track_point_t *)b;
   return g_date_time_compare(pa->time, pb->time);
 }
 
-static gint _sort_segment(gconstpointer a, gconstpointer b)
+static gint _sort_segment(const gconstpointer a, const gconstpointer b)
 {
   const dt_gpx_track_segment_t *pa = (const dt_gpx_track_segment_t *)a;
   const dt_gpx_track_segment_t *pb = (const dt_gpx_track_segment_t *)b;
@@ -150,7 +150,7 @@ void dt_gpx_destroy(dt_gpx_t *gpx)
   g_free(gpx);
 }
 
-gboolean dt_gpx_get_location(dt_gpx_t *gpx, GDateTime *timestamp, dt_image_geoloc_t *geoloc)
+gboolean dt_gpx_get_location(const dt_gpx_t *gpx, GDateTime *timestamp, dt_image_geoloc_t *geoloc)
 {
   g_assert(gpx != NULL);
 
@@ -159,7 +159,7 @@ gboolean dt_gpx_get_location(dt_gpx_t *gpx, GDateTime *timestamp, dt_image_geolo
 
   for(GList *item = gpx->trkpts; item; item = g_list_next(item))
   {
-    dt_gpx_track_point_t *tp = item->data;
+    const dt_gpx_track_point_t *tp = item->data;
 
     /* if timestamp is out of time range return false but fill
        closest location value start or end point */
@@ -172,13 +172,13 @@ gboolean dt_gpx_get_location(dt_gpx_t *gpx, GDateTime *timestamp, dt_image_geolo
       return FALSE;
     }
 
-    dt_gpx_track_point_t *tp_next = item->next->data;
+    const dt_gpx_track_point_t *tp_next = item->next->data;
     /* check if timestamp is within current and next trackpoint */
     const gint cmp_n = g_date_time_compare(timestamp, tp_next->time);
     if(item->next && cmp_n <= 0)
     {
-      GTimeSpan seg_diff = g_date_time_difference(tp_next->time, tp->time);
-      GTimeSpan diff = g_date_time_difference(timestamp, tp->time);
+      const GTimeSpan seg_diff = g_date_time_difference(tp_next->time, tp->time);
+      const GTimeSpan diff = g_date_time_difference(timestamp, tp->time);
       if(seg_diff == 0 || diff == 0)
       {
         geoloc->longitude = tp->longitude;
@@ -264,7 +264,7 @@ static void _gpx_parse_error(GError **error)
  */
 void _gpx_parser_start_element(GMarkupParseContext *ctx, const gchar *element_name,
                                const gchar **attribute_names, const gchar **attribute_values,
-                               gpointer user_data, GError **error)
+                               const gpointer user_data, GError **error)
 {
   g_return_if_fail(*error == NULL);
 
@@ -374,7 +374,7 @@ void _gpx_parser_start_element(GMarkupParseContext *ctx, const gchar *element_na
   _gpx_parse_error(error);
 }
 
-void _gpx_parser_end_element(GMarkupParseContext *context, const gchar *element_name, gpointer user_data,
+void _gpx_parser_end_element(GMarkupParseContext *context, const gchar *element_name, const gpointer user_data,
                              GError **error)
 {
   g_return_if_fail(*error == NULL);
@@ -407,7 +407,7 @@ void _gpx_parser_end_element(GMarkupParseContext *context, const gchar *element_
   }
 }
 
-void _gpx_parser_text(GMarkupParseContext *context, const gchar *text, gsize text_len, gpointer user_data,
+void _gpx_parser_text(GMarkupParseContext *context, const gchar *text, gsize text_len, const gpointer user_data,
                       GError **error)
 {
   g_return_if_fail(*error == NULL);
@@ -458,26 +458,26 @@ void _gpx_parser_text(GMarkupParseContext *context, const gchar *text, gsize tex
     gpx->current_track_point->elevation = g_ascii_strtod(text, NULL);
 }
 
-GList *dt_gpx_get_trkseg(struct dt_gpx_t *gpx)
+GList *dt_gpx_get_trkseg(const struct dt_gpx_t *gpx)
 {
   return (gpx != NULL)? gpx->trksegs
                       : NULL;
 }
 
-GList *dt_gpx_get_trkpts(dt_gpx_t *gpx, const guint segid)
+GList *dt_gpx_get_trkpts(const dt_gpx_t *gpx, const guint segid)
 {
   if(gpx == NULL)
     return NULL;
 
   GList *pts = NULL;
   GList *ts = g_list_nth(gpx->trksegs, segid);
   if(!ts) return pts;
-  dt_gpx_track_segment_t *tsd = ts->data;
+  const dt_gpx_track_segment_t *tsd = ts->data;
   GList *tps = g_list_find(gpx->trkpts, tsd->trkpt);
   if(!tps) return pts;
   for(GList *tp = tps; tp; tp = g_list_next(tp))
   {
-    dt_gpx_track_point_t *tpd = tp->data;
+    const dt_gpx_track_point_t *tpd = tp->data;
     if(tpd->segid != segid) return pts;
     dt_geo_map_display_point_t *p = g_malloc0(sizeof(dt_geo_map_display_point_t));
     p->lat = tpd->latitude;
@@ -491,14 +491,14 @@ GList *dt_gpx_get_trkpts(dt_gpx_t *gpx, const guint segid)
  * Geodesic interpolation functions
  * ------------------------------------------------------------------------*/
 
-void dt_gpx_geodesic_distance(double lat1, double lon1,
-                              double lat2, double lon2,
+void dt_gpx_geodesic_distance(const double lat1, const double lon1,
+                              const double lat2, const double lon2,
                               double *d, double *delta)
 {
-  const double lat_rad_1 = lat1 * M_PI / 180;
-  const double lat_rad_2 = lat2 * M_PI / 180;
-  const double lon_rad_1 = lon1 * M_PI / 180;
-  const double lon_rad_2 = lon2 * M_PI / 180;
+  const double lat_rad_1 = deg2rad(lat1);
+  const double lat_rad_2 = deg2rad(lat2);
+  const double lon_rad_1 = deg2rad(lon1);
+  const double lon_rad_2 = deg2rad(lon2);
   const double delta_lat_rad = lat_rad_2 - lat_rad_1;
   const double delta_lon_rad = lon_rad_2 - lon_rad_1;
   const double sin_delta_lat_rad = sin(delta_lat_rad / 2);
@@ -516,7 +516,7 @@ void dt_gpx_geodesic_intermediate_point(const double lat1, const double lon1,
                                         const double lat2, const double lon2,
                                         const double delta,
                                         const gboolean first_time,
-                                        double f,
+                                        const double f,
                                         double *lat, double *lon)
 {
   static double lat_rad_1;
@@ -535,16 +535,16 @@ void dt_gpx_geodesic_intermediate_point(const double lat1, const double lon1,
 
   if(first_time)
   {
-    lat_rad_1 = lat1 * M_PI / 180;
+    lat_rad_1 = deg2rad(lat1);
     sin_lat_rad_1 = sin(lat_rad_1);
     cos_lat_rad_1 = cos(lat_rad_1);
-    lat_rad_2 = lat2 * M_PI / 180;
+    lat_rad_2 = deg2rad(lat2);
     sin_lat_rad_2 = sin(lat_rad_2);
     cos_lat_rad_2 = cos(lat_rad_2);
-    lon_rad_1 = lon1 * M_PI / 180;
+    lon_rad_1 = deg2rad(lon1);
     sin_lon_rad_1 = sin(lon_rad_1);
     cos_lon_rad_1 = cos(lon_rad_1);
-    lon_rad_2 = lon2 * M_PI / 180;
+    lon_rad_2 = deg2rad(lon2);
     sin_lon_rad_2 = sin(lon_rad_2);
     cos_lon_rad_2 = cos(lon_rad_2);
     sin_delta = sin(delta);
@@ -558,8 +558,8 @@ void dt_gpx_geodesic_intermediate_point(const double lat1, const double lon1,
   const double lat_rad = atan2(z, sqrt(x * x + y * y)); /* latitude of intermediate point in radians */
   const double lon_rad = atan2(y, x);                   /* longitude of intermediate point in radians */
 
-  *lat = lat_rad / M_PI * 180;
-  *lon = lon_rad / M_PI * 180;
+  *lat = rad2deg(lat_rad);
+  *lon = rad2deg(lon_rad);
 }
 /* -------- end of Geodesic interpolation functions -----------------------*/
 

src/common/gpx.h

@@ -53,13 +53,13 @@ void dt_gpx_destroy(struct dt_gpx_t *gpx);
 /* fetch the lon,lat coords for time t, if within time range
   of gpx record return TRUE, FALSE is returned if out of time frame
   and closest record of lon,lat is filled */
-gboolean dt_gpx_get_location(struct dt_gpx_t *, GDateTime *timestamp, dt_image_geoloc_t *geoloc);
+gboolean dt_gpx_get_location(const struct dt_gpx_t *, GDateTime *timestamp, dt_image_geoloc_t *geoloc);
 
 // get the list of track segments
-GList *dt_gpx_get_trkseg(struct dt_gpx_t *gpx);
+GList *dt_gpx_get_trkseg(const struct dt_gpx_t *gpx);
 
 // get the list of track points for a track segment
-GList *dt_gpx_get_trkpts(struct dt_gpx_t *gpx, const guint segid);
+GList *dt_gpx_get_trkpts(const struct dt_gpx_t *gpx, const guint segid);
 
 // get the distance on the geodesic line and the angle delta referred to earth center
 void dt_gpx_geodesic_distance(double lat1, double lon1,