star-cem

Compute the position of the sun
git clone git://git.meso-star.com/star-cem.git
Log | Files | Refs | README | LICENSE

commit 23fdcc4c2f8feb5c2a1e43c3fa6aa95db904fbd7
parent 07eb15d702237c210ef3cd90092dbfc8ae0e9e84
Author: Christophe Coustet <christophe.coustet@meso-star.com>
Date:   Thu, 18 Jun 2026 13:42:50 +0200

Add an utility to convert sun angles to sun vector

Also improve wording an the API comments and the man page.

Diffstat:
Mdoc/scem.1 | 75++++++++++++++++++++++++++++++++++++++++++++++-----------------------------
Msrc/scem.c | 30++++++++++++++++++++++++++++++
Msrc/scem.h | 16++++++++++++----
Msrc/scem_main.c | 27++++++++++++++++++++-------
Msrc/test_scem_sun_position.c | 7+++++++
5 files changed, 115 insertions(+), 40 deletions(-)

diff --git a/doc/scem.1 b/doc/scem.1 @@ -23,7 +23,7 @@ .\"""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""" .Sh SYNOPSIS .Nm -.Op Fl hvr +.Op Fl DhVvr .Op Fl a Ar algo .Op Fl d Ar utc_date .Ar latitude @@ -34,36 +34,43 @@ computes the position of the sun relative to a given point on the Earth's surface, on a given date. .Pp -The latitude must be in -.Bq -90, 90 -degrees relative to the equator. +The latitude must be in [-90, 90] degrees relative to the equator. It is positive towards the north. -The longitude must be in -.Bq -180, 180 -degrees relative to Greenwich. +The longitude must be in [-180, 180] degrees relative to Greenwich. The angle is positive towards the east. .Pp -On output, -.Nm -displays the position of the sun in elevation and azimuth on the standard -output. -The message is formatted as follows: +The output depends on two options, +.Fl D +and +.Fl V . +If +.Fl D +is not provided, +.No +outputs the position of the sun in elevation and azimuth to the standard +output with the following format: .Bd -literal -offset Ds "%f %f\\n", elevation, azimuth .Ed .Pp -The elevation angle is in -.Bq -90, 90 -degrees. +The elevation angle is in [-90, 90] degrees. It defines the elevation of the sun relative to the horizon. -The azimuth angle is in -.Bq 0, 360 -degrees, with 0° indicating north, +The azimuth angle is in [0, 360] degrees, with 0° indicating north, 90° east, 180° south, and 270° west. .Pp -The options are as follows +If +.Fl V +is provided, +.Nm +outputs the direction from the sun as a 3D vector to the standard output +with the following format: +.Bd -literal -offset Ds +"%f %f %f\\n", vec[0], vec[1], vec[2] +.Ed +.Pp +The options are as follows: .Bl -tag -width Ds .\"""""""""""""""""""""""""""""""""""""" .It Fl a Ar algo @@ -75,13 +82,15 @@ algorithm. The solar position algorithms are as follows: .Bl -tag -width Ds .It Cm meeus -The algorithm described in the Jean Meeus' book, -.Dq Astronomical Algorithm . +The algorithm described in Jean Meeus' book, +.Dq Astronomical Algorithms . .It Cm psa The algorithm developed by the .Dq Plataforma Solar de Almería . .El .\"""""""""""""""""""""""""""""""""""""" +.It Fl D +Suppress the sun direction from the output. .It Fl d Ar utc_date Date on which the position of the sun is calculated at the specified location. @@ -99,29 +108,37 @@ date -u +"%Y-%m-%dT%H:%M:%S" .It Fl h Output short help and exit. .\"""""""""""""""""""""""""""""""""""""" +.It Fl V +Add the 3D vector from the sun to the observer to the output with the first, +second and third components pointing east, north and up, respectively. .It Fl v Make .Nm verbose. .\"""""""""""""""""""""""""""""""""""""" .It Fl r -The solar output coordinates are defined in radians rather than degrees. +Output the solar direction in radians rather than degrees. +Meaningless if +.Fl D +is defined. .El .\"""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""" -.Sh EXIT +.Sh EXIT STATUS .Ex -std .\"""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""" .Sh EXAMPLES -Calculate the position of the sun on the current date for a latitude of -43.559962° east and a longitude of 1.468150° north: +Compute the current position of the sun for a latitude of 43.559962° north and a +longitude of 1.468150° east. +Output both the direction in radians and the corresponding 3D vector: .Bd -literal -offset Ds -scem 43.559962 1.468150 +scem -V -r 43.559962 1.468150 .Ed .Pp -Same as above, but setting the observation date to October 21, 2015, at -7:28:00 a.m.: +Same location as above, setting the observation date to October 21, 2015, at +11:28:00 AM UTC+2. +Output the 3D vector only: .Bd -literal -offset Ds -scem -d 2015-10-21T7:28:00 43.559962 1.468150 +scem -D -V -d 2015-10-21T09:28:00 43.559962 1.468150 .Ed .\"""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""" .Sh SEE ALSO diff --git a/src/scem.c b/src/scem.c @@ -16,6 +16,8 @@ #include "scem.h" #include "scem_c.h" +#include <math.h> + /******************************************************************************* * Exported functions ******************************************************************************/ @@ -48,3 +50,31 @@ scem_sun_position_from_earth } return res; } + +res_T +scem_sun_position_to_sun_vector + (const struct scem_sun_pos* spherical, + double sun_dir[3]) +{ + double cos_azimuth; + double sin_azimuth; + double cos_elevation; + double sin_elevation; + double azimuth; + double elevation; + + if(!spherical || !sun_dir) return RES_BAD_ARG; + + azimuth = spherical->azimuth; + elevation = spherical->elevation; + + cos_azimuth = cos(azimuth); + sin_azimuth = sin(azimuth); + cos_elevation = cos(elevation); + sin_elevation = sin(elevation); + + sun_dir[0] = -(cos_elevation * sin_azimuth); + sun_dir[1] = -(cos_elevation * cos_azimuth); + sun_dir[2] = -(sin_elevation); + return RES_OK; +} diff --git a/src/scem.h b/src/scem.h @@ -73,10 +73,10 @@ scem_sun_algo_cstr(const enum scem_sun_algo algo) BEGIN_DECLS -/* Compute the position of the Sun in local coordinates (elevation and azimuth) - * in radians for a given observer at the surface of the Earth (can be extended - * to other * planets if required) i.e. the observer has to be located at the - * surface for the required date */ +/* For a given location and time at the surface of the Earth (the observer has + * to be located at the surface), compute the position of the sun in local + * coordinates (elevation and azimuth) in radians. + * Could be extended to other planets if required. */ SCEM_API res_T scem_sun_position_from_earth (const struct tm* time, /* In UTC */ @@ -84,6 +84,14 @@ scem_sun_position_from_earth const enum scem_sun_algo algo, struct scem_sun_pos* sun); +/* Convert the position of the sun (elevation and azimuth) in radians, as + * produced by scem_sun_position_from_earth, to a 3D vector (sun to observer). + * First component east, second component north, third component up. */ +SCEM_API res_T +scem_sun_position_to_sun_vector + (const struct scem_sun_pos* sun, + double sun_vector[3]); + END_DECLS #endif /* SCEM_H */ diff --git a/src/scem_main.c b/src/scem_main.c @@ -18,6 +18,7 @@ #include "scem.h" +#include <rsys/rsys.h> #include <rsys/cstr.h> #include <rsys/mem_allocator.h> @@ -33,9 +34,11 @@ struct args { int radian; /* Sun position in radians instead of degrees */ int verbose; int quit; + int no_dir_output; + int vec_output; }; static const struct args ARGS_DEFAULT = { - SCEM_LOCATION_NULL__, {0}, SCEM_SUN_MEEUS, 0, 0, 0 + SCEM_LOCATION_NULL__, {0}, SCEM_SUN_MEEUS, 0, 0, 0, 0, 0 }; /******************************************************************************* @@ -45,7 +48,7 @@ static INLINE void usage(FILE* stream) { fprintf(stream, - "usage: scem [-hrv] [-a algo] [-d utc_date] latitude longitude\n"); + "usage: scem [-DhrVv] [-a algo] [-d utc_date] latitude longitude\n"); } static res_T @@ -125,11 +128,13 @@ args_init(struct args* args, int argc, char** argv) if((res = get_current_date(&args->time)) != RES_OK) goto error; /* Don't process the last 2 args as options, they are lat and long */ - while((opt = getopt(argc - 2, argv, "a:d:hvr")) != -1) { + while((opt = getopt(argc - 2, argv, "a:Dd:hVvr")) != -1) { switch(opt) { case 'a': res = parse_algo(optarg, &args->algo); break; + case 'D': args->no_dir_output = 1; break; case 'd': res = parse_date(optarg, &args->time); break; case 'h': usage(stdout); args->quit = 1; goto exit; + case 'V': args->vec_output = 1; break; case 'v': args->verbose = 1; break; case 'r': args->radian = 1; break; default: res = RES_BAD_ARG; break; @@ -165,6 +170,7 @@ main(int argc, char** argv) { struct args args = ARGS_DEFAULT; struct scem_sun_pos pos = SCEM_SUN_POS_NULL; + double sun_dir[3]; int err = 0; res_T res = RES_OK; @@ -180,10 +186,17 @@ main(int argc, char** argv) goto error; } - if(args.radian) { - printf("%g %g\n", pos.elevation, pos.azimuth); - } else { - printf("%g %g\n", MRAD2DEG(pos.elevation), MRAD2DEG(pos.azimuth)); + if(!args.no_dir_output) { + if(args.radian) { + printf("%g %g\n", pos.elevation, pos.azimuth); + } else { + printf("%g %g\n", MRAD2DEG(pos.elevation), MRAD2DEG(pos.azimuth)); + } + } + if(args.vec_output) { + res = scem_sun_position_to_sun_vector(&pos, sun_dir); + if(res != RES_OK) goto error; + fprintf(stderr, "%g %g %g\n", SPLIT3(sun_dir)); } exit: diff --git a/src/test_scem_sun_position.c b/src/test_scem_sun_position.c @@ -33,6 +33,7 @@ test_api(void) struct tm utc = {0}; struct scem_location pos = SCEM_LOCATION_NULL; struct scem_sun_pos sun = SCEM_SUN_POS_NULL; + double vec[3]; utc.tm_sec = 0; utc.tm_min = 28; @@ -47,6 +48,10 @@ test_api(void) pos.longitude = 1.468150; printf("Lat: %g°; lon: %g°\n", pos.latitude, pos.longitude); + CHK(scem_sun_position_to_sun_vector(NULL, NULL) == RES_BAD_ARG); + CHK(scem_sun_position_to_sun_vector(&sun, NULL) == RES_BAD_ARG); + CHK(scem_sun_position_to_sun_vector(NULL, vec) == RES_BAD_ARG); + CHK(scem_sun_position_from_earth(NULL, &pos, algo, &sun) == RES_BAD_ARG); CHK(scem_sun_position_from_earth(&utc, NULL, algo, &sun) == RES_BAD_ARG); CHK(scem_sun_position_from_earth(&utc, &pos, SCEM_SUN_ALGO_NONE__, &sun) @@ -71,6 +76,8 @@ test_api(void) scem_sun_algo_cstr(algo), MRAD2DEG(sun.elevation), MRAD2DEG(sun.azimuth)); + CHK(scem_sun_position_to_sun_vector(&sun, vec) == RES_OK); + printf("Vector: %g %g %g\n", SPLIT3(vec)); } }