1 /*
2 swetest.c A test program
3
4 Authors: Dieter Koch and Alois Treindl, Astrodienst Zuerich
5
6 **************************************************************/
7
8 /* Copyright (C) 1997 - 2021 Astrodienst AG, Switzerland. All rights reserved.
9
10 License conditions
11 ------------------
12
13 This file is part of Swiss Ephemeris.
14
15 Swiss Ephemeris is distributed with NO WARRANTY OF ANY KIND. No author
16 or distributor accepts any responsibility for the consequences of using it,
17 or for whether it serves any particular purpose or works at all, unless he
18 or she says so in writing.
19
20 Swiss Ephemeris is made available by its authors under a dual licensing
21 system. The software developer, who uses any part of Swiss Ephemeris
22 in his or her software, must choose between one of the two license models,
23 which are
24 a) GNU Affero General Public License (AGPL)
25 b) Swiss Ephemeris Professional License
26
27 The choice must be made before the software developer distributes software
28 containing parts of Swiss Ephemeris to others, and before any public
29 service using the developed software is activated.
30
31 If the developer choses the AGPL software license, he or she must fulfill
32 the conditions of that license, which includes the obligation to place his
33 or her whole software project under the AGPL or a compatible license.
34 See https://www.gnu.org/licenses/agpl-3.0.html
35
36 If the developer choses the Swiss Ephemeris Professional license,
37 he must follow the instructions as found in http://www.astro.com/swisseph/
38 and purchase the Swiss Ephemeris Professional Edition from Astrodienst
39 and sign the corresponding license contract.
40
41 The License grants you the right to use, copy, modify and redistribute
42 Swiss Ephemeris, but only under certain conditions described in the License.
43 Among other things, the License requires that the copyright notices and
44 this notice be preserved on all copies.
45
46 Authors of the Swiss Ephemeris: Dieter Koch and Alois Treindl
47
48 The authors of Swiss Ephemeris have no control or influence over any of
49 the derived works, i.e. over software or services created by other
50 programmers which use Swiss Ephemeris functions.
51
52 The names of the authors or of the copyright holder (Astrodienst) must not
53 be used for promoting any software, product or service which uses or contains
54 the Swiss Ephemeris. This copyright notice is the ONLY place where the
55 names of the authors can legally appear, except in cases where they have
56 given special permission in writing.
57
58 The trademarks 'Swiss Ephemeris' and 'Swiss Ephemeris inside' may be used
59 for promoting such software, products or services.
60 */
61
62 /* attention: Microsoft Compiler does not accept strings > 2048 char */
63
64 static char *infocmd0 = "\n\
65 Swetest computes a complete set of geocentric planetary positions,\n\
66 for a given date or a sequence of dates.\n\
67 Input can either be a date or an absolute julian day number.\n\
68 0:00 (midnight).\n\
69 With the proper options, swetest can be used to output a printed\n\
70 ephemeris and transfer the data into other programs like spreadsheets\n\
71 for graphical display.\n\
72 Version: \n\
73 \n";
74 static char *infocmd1 = "\n\
75 Command line options:\n\
76 help commands:\n\
77 -?, -h display whole info\n\
78 -hcmd display commands\n\
79 -hplan display planet numbers\n\
80 -hform display format characters\n\
81 -hdate display input date format\n\
82 -hexamp display examples\n\
83 -glp report file location of library\n\
84 input time formats:\n\
85 -bDATE begin date; e.g. -b1.1.1992 if\n\
86 Note: the date format is day month year (European style).\n\
87 -bj... begin date as an absolute Julian day number; e.g. -bj2415020.5\n\
88 -j... same as -bj\n\
89 -tHH.MMSS input time (as Ephemeris Time)\n\
90 -ut input date is Universal Time (UT1)\n\
91 -utHH:MM:SS input time (as Universal Time)\n\
92 -utHH.MMSS input time (as Universal Time)\n\
93 -utcHH.MM:SS input time (as Universal Time Coordinated UTC)\n\
94 output time for eclipses, occultations, risings/settings is UT by default\n\
95 -lmt output date/time is LMT (with -geopos)\n\
96 -lat output date/time is LAT (with -geopos)\n\
97 object, number of steps, step with\n\
98 -pSEQ planet sequence to be computed.\n\
99 See the letter coding below.\n\
100 -dX differential ephemeris: print differential ephemeris between\n\
101 body X and each body in list given by -p\n\
102 example: -p2 -d0 -fJl -n366 -b1.1.1992 prints the longitude\n\
103 distance between SUN (planet 0) and MERCURY (planet 2)\n\
104 for a full year starting at 1 Jan 1992.\n\
105 -dhX differential ephemeris: print differential ephemeris between\n\
106 heliocentric body X and each body in list given by -p\n\
107 example: -p8 -dh8 -ftl -n36600 -b1.1.1500 -s5 prints the longitude\n\
108 distance between geocentric and heliocentric Neptune (planet 8)\n\
109 for 500 year starting at 1 Jan 1500.\n\
110 Using this option mostly makes sense for a single planet\n\
111 to find out how much its geocentric and heliocentric positions can differ\n\
112 over extended periods of time\n\
113 -DX midpoint ephemeris, works the same way as the differential\n\
114 mode -d described above, but outputs the midpoint position.\n\
115 -nN output data for N consecutive timesteps; if no -n option\n\
116 is given, the default is 1. If the option -n without a\n\
117 number is given, the default is 20.\n\
118 -sN timestep N days, default 1. This option is only meaningful\n\
119 when combined with option -n.\n\
120 If an 'y' is appended, the time step is in years instead of days, \n\
121 for example -s10y for a time step of 10 years.\n\
122 If an 'mo' is appended, the time step is in months instead of days, \n\
123 for example -s3mo for a time step of 3 months.\n\
124 If an 'm' is appended, the time step is in minutes instead of days, \n\
125 for example -s15m for a time step of 15 minutes.\n\
126 If an 's' is appended, the time step is in seconds instead of days, \n\
127 for example -s1s for a time step of 1 second.\n\
128 ";
129 static char *infocmd2 = "\
130 output format:\n\
131 -fSEQ use SEQ as format sequence for the output columns;\n\
132 default is PLBRS.\n\
133 -head don\'t print the header before the planet data. This option\n\
134 is useful when you want to paste the output into a\n\
135 spreadsheet for displaying graphical ephemeris.\n\
136 +head header before every step (with -s..) \n\
137 -gPPP use PPP as gap between output columns; default is a single\n\
138 blank. -g followed by white space sets the\n\
139 gap to the TAB character; which is useful for data entry\n\
140 into spreadsheets.\n\
141 -hor list data for multiple planets 'horizontally' in same line.\n\
142 all columns of -fSEQ are repeated except time colums tTJyY.\n\
143 astrological house system:\n\
144 -house[long,lat,hsys] \n\
145 include house cusps. The longitude, latitude (degrees with\n\
146 DECIMAL fraction) and house system letter can be given, with\n\
147 commas separated, + for east and north. If none are given,\n\
148 Greenwich UK and Placidus is used: 0.00,51.50,p.\n\
149 The output lists 12 house cusps, Asc, MC, ARMC, Vertex,\n\
150 Equatorial Ascendant, co-Ascendant as defined by Walter Koch, \n\
151 co-Ascendant as defined by Michael Munkasey, and Polar Ascendant. \n\
152 Houses can only be computed if option -ut is given.\n\
153 A equal\n\
154 B Alcabitius\n\
155 C Campanus\n\
156 D equal / MC\n\
157 E equal = A\n\
158 F Carter poli-equatorial\n\
159 G 36 Gauquelin sectors\n\
160 H horizon / azimuth\n\
161 I Sunshine\n\
162 i Sunshine alternative\n\
163 K Koch\n\
164 L Pullen S-delta\n\
165 M Morinus\n\
166 N Whole sign, Aries = 1st house\n\
167 O Porphyry\n\
168 P Placidus\n\
169 Q Pullen S-ratio\n\
170 R Regiomontanus\n\
171 S Sripati\n\
172 T Polich/Page (\"topocentric\")\n\
173 U Krusinski-Pisa-Goelzer\n\
174 V equal Vehlow\n\
175 W equal, whole sign\n\
176 X axial rotation system/ Meridian houses\n\
177 Y APC houses\n\
178 The use of lower case letters is deprecated. They will have a\n\
179 different meaning in future releases of Swiss Ephemeris.\n\
180 -hsy[hsys] \n\
181 house system to be used (for house positions of planets)\n\
182 for long, lat, hsys, see -house\n\
183 The use of lower case letters is deprecated. They will have a\n\
184 different meaning in future releases of Swiss Ephemeris.\n\
185 ";
186 static char *infocmd3 = "\
187 -geopos[long,lat,elev] \n\
188 Geographic position. Can be used for azimuth and altitude\n\
189 or house cusps calculations.\n\
190 The longitude, latitude (degrees with DECIMAL fraction)\n\
191 and elevation (meters) can be given, with\n\
192 commas separated, + for east and north. If none are given,\n\
193 Greenwich is used: 0,51.5,0.\n\
194 For topocentric planet positions please user the parameter -topo\n\
195 sidereal astrology:\n\
196 -ay.. ayanamsha, with number of method, e.g. ay0 for Fagan/Bradley\n\
197 -sid.. sidereal, with number of method (see below)\n\
198 -sidt0.. dito, but planets are projected on the ecliptic plane of the\n\
199 reference date of the ayanamsha (more info in general documentation\n\
200 www.astro.com/swisseph/swisseph.htm)\n\
201 -sidsp.. dito, but planets are projected on the solar system plane.\n\
202 (see www.astro.com/swisseph/swisseph.htm)\n\
203 -sidudef[jd,ay0,...] sidereal, with user defined ayanamsha; \n\
204 jd=julian day number in TT/ET\n\
205 ay0=initial value of ayanamsha, \n\
206 ...=optional parameters, comma-sparated:\n\
207 'jdisut': ayanamsha reference date is UT\n\
208 'eclt0': project on ecliptic of reference date (like -sidt0..)\n\
209 'ssyplane': project on solar system plane (like -sidsp..)\n\
210 e.g. '-sidudef2452163.8333333,25.0,jdisut': ayanamsha is 25.0° on JD 2452163.8333333 UT\n\
211 number of ayanamsha method:\n\
212 0 for Fagan/Bradley\n\
213 1 for Lahiri\n\
214 2 for De Luce\n\
215 3 for Raman\n\
216 4 for Usha/Shashi\n\
217 5 for Krishnamurti\n\
218 6 for Djwhal Khul\n\
219 7 for Yukteshwar\n\
220 8 for J.N. Bhasin\n\
221 9 for Babylonian/Kugler 1\n\
222 10 for Babylonian/Kugler 2\n\
223 11 for Babylonian/Kugler 3\n\
224 12 for Babylonian/Huber\n\
225 13 for Babylonian/Eta Piscium\n\
226 14 for Babylonian/Aldebaran = 15 Tau\n\
227 15 for Hipparchos\n\
228 16 for Sassanian\n\
229 17 for Galact. Center = 0 Sag\n\
230 18 for J2000\n\
231 19 for J1900\n\
232 20 for B1950\n\
233 21 for Suryasiddhanta\n\
234 22 for Suryasiddhanta, mean Sun\n\
235 23 for Aryabhata\n\
236 24 for Aryabhata, mean Sun\n\
237 25 for SS Revati\n\
238 26 for SS Citra\n\
239 27 for True Citra\n\
240 28 for True Revati\n\
241 29 for True Pushya (PVRN Rao)\n\
242 30 for Galactic (Gil Brand)\n\
243 31 for Galactic Equator (IAU1958)\n\
244 32 for Galactic Equator\n\
245 33 for Galactic Equator mid-Mula\n\
246 34 for Skydram (Mardyks)\n\
247 35 for True Mula (Chandra Hari)\n\
248 36 Dhruva/Gal.Center/Mula (Wilhelm)\n\
249 37 Aryabhata 522\n\
250 38 Babylonian/Britton\n\
251 39 Vedic/Sheoran\n\
252 40 Cochrane (Gal.Center = 0 Cap)\n\
253 41 Galactic Equator (Fiorenza)\n\
254 42 Vettius Valens\n\
255 43 Lahiri 1940\n\
256 44 Lahiri VP285 (1980)\n\
257 45 Krishnamurti VP291\n\
258 46 Lahiri ICRC\n\
259 ephemeris specifications:\n\
260 -edirPATH change the directory of the ephemeris files \n\
261 -eswe swiss ephemeris\n\
262 -ejpl jpl ephemeris (DE431), or with ephemeris file name\n\
263 -ejplde200.eph \n\
264 -emos moshier ephemeris\n\
265 -true true positions\n\
266 -noaberr no aberration\n\
267 -nodefl no gravitational light deflection\n\
268 -noaberr -nodefl astrometric positions\n\
269 -j2000 no precession (i.e. J2000 positions)\n\
270 -icrs ICRS (use Internat. Celestial Reference System)\n\
271 -nonut no nutation \n\
272 ";
273 static char *infocmd4 = "\
274 -speed calculate high precision speed \n\
275 -speed3 'low' precision speed from 3 positions \n\
276 do not use this option. -speed parameter\n\
277 is faster and more precise \n\
278 -iXX force iflag to value XX\n\
279 -testaa96 test example in AA 96, B37,\n\
280 i.e. venus, j2450442.5, DE200.\n\
281 attention: use precession IAU1976\n\
282 and nutation 1980 (s. swephlib.h)\n\
283 -testaa95\n\
284 -testaa97\n\
285 -roundsec round to seconds\n\
286 -roundmin round to minutes\n\
287 -ep use extra precision in output for some data\n\
288 -dms use dms instead of fractions, at some places\n\
289 -lim print ephemeris file range\n\
290 observer position:\n\
291 -hel compute heliocentric positions\n\
292 -bary compute barycentric positions (bar. earth instead of node) \n\
293 -topo[long,lat,elev] \n\
294 topocentric positions. The longitude, latitude (degrees with\n\
295 DECIMAL fraction) and elevation (meters) can be given, with\n\
296 commas separated, + for east and north. If none are given,\n\
297 Greenwich is used 0.00,51.50,0\n\
298 -pc... compute planetocentric positions\n\
299 to specify the central body, use the internal object number\n\
300 of Swiss Ephemeris, e.g. 3 for Venus, 4 for Mars, \n\
301 -pc3 Venus-centric \n\
302 -pc4 Mars-centric \n\
303 -pc5 Jupiter-centric (barycenter)\n\
304 -pc9599 Jupiter-centric (center of body)\n\
305 -pc9699 Saturn-centric (center of body)\n\
306 For asteroids use MPC number + 10000, e.g.\n\
307 -pc10433 Eros-centric (Eros = 433 + 10000)\n\
308 orbital elements:\n\
309 -orbel compute osculating orbital elements relative to the\n\
310 mean ecliptic J2000. (Note, all values, including time of\n\
311 pericenter vary considerably depending on the date for which the\n\
312 osculating ellipse is calculated\n\
313 \n\
314 special events:\n\
315 -solecl solar eclipse\n\
316 output 1st line:\n\
317 eclipse date,\n\
318 time of maximum (UT):\n\
319 geocentric angle between centre of Sun and Moon reaches minimum.\n\
320 core shadow width (negative with total eclipses),\n\
321 eclipse magnitudes:\n\
322 1. NASA method (= 2. with partial ecl. and \n\
323 ratio lunar/solar diameter with total and annular ecl.)\n\
324 2. fraction of solar diameter covered by moon;\n\
325 if the value is > 1, it means that Moon covers more than\n\
326 just the solar disk\n\
327 3. fraction of solar disc covered by moon (obscuration)\n\
328 with total and annular eclipses it is the ratio of\n\
329 the sizes of the solar disk and the lunar disk.\n\
330 Saros series and eclipse number\n\
331 Julian day number (6-digit fraction) of maximum\n\
332 output 2nd line:\n\
333 start and end times for partial and total phases\n\
334 delta t in sec\n\
335 output 3rd line:\n\
336 geographical longitude and latitude of maximum eclipse,\n\
337 totality duration at that geographical position,\n\
338 output with -local, see below.\n\
339 -occult occultation of planet or star by the moon. Use -p to \n\
340 specify planet (-pf -xfAldebaran for stars) \n\
341 output format same as with -solecl, with the following differences:\n\
342 Magnitude is defined like no. 2. with solar eclipses.\n\
343 There are no saros series.\n\
344 ";
345 static char *infocmd5 = "\
346 -lunecl lunar eclipse\n\
347 output 1st line:\n\
348 eclipse date,\n\
349 time of maximum (UT),\n\
350 eclipse magnitudes: umbral and penumbral\n\
351 method as method 2 with solar eclipses\n\
352 Saros series and eclipse number \n\
353 Julian day number (6-digit fraction) of maximum\n\
354 output 2nd line:\n\
355 6 contacts for start and end of penumbral, partial, and\n\
356 total phase\n\
357 delta t in sec\n\
358 output 3rd line:\n\
359 geographic position where the Moon is in zenith at maximum eclipse\n\
360 -local only with -solecl or -occult, if the next event of this\n\
361 kind is wanted for a given geogr. position.\n\
362 Use -geopos[long,lat,elev] to specify that position.\n\
363 If -local is not set, the program \n\
364 searches for the next event anywhere on earth.\n\
365 output 1st line:\n\
366 eclipse date,\n\
367 time of maximum,\n\
368 eclipse magnitudes, as with global solar eclipse function \n\
369 (with occultations: only diameter method, see solar eclipses, method 2)\n\
370 Saros series and eclipse number (with solar eclipses only)\n\
371 Julian day number (6-digit fraction) of maximum\n\
372 output 2nd line:\n\
373 local eclipse duration for totality (zero with partial occultations)\n\
374 local four contacts,\n\
375 delta t in sec\n\
376 Occultations with the remark \"(daytime)\" cannot be observed because\n\
377 they are taking place by daylight. Occultations with the remark\n\
378 \"(sunrise)\" or \"(sunset)\" can be observed only partly because part\n\
379 of them takes place in daylight.\n\
380 -hev[type] heliacal events,\n\
381 type 1 = heliacal rising\n\
382 type 2 = heliacal setting\n\
383 type 3 = evening first\n\
384 type 4 = morning last\n\
385 type 0 or missing = all four events are listed.\n\
386 -rise rising and setting of a planet or star.\n\
387 Use -geopos[long,lat,elev] to specify geographical position.\n\
388 -metr southern and northern meridian transit of a planet of star\n\
389 Use -geopos[long,lat,elev] to specify geographical position.\n\
390 specifications for eclipses:\n\
391 -total total eclipse (only with -solecl, -lunecl)\n\
392 -partial partial eclipse (only with -solecl, -lunecl)\n\
393 -annular annular eclipse (only with -solecl)\n\
394 -anntot annular-total (hybrid) eclipse (only with -solecl)\n\
395 -penumbral penumbral lunar eclipse (only with -lunecl)\n\
396 -central central eclipse (only with -solecl, nonlocal)\n\
397 -noncentral non-central eclipse (only with -solecl, nonlocal)\n\
398 ";
399 static char *infocmd6 = "\
400 specifications for risings and settings:\n\
401 -norefrac neglect refraction (with option -rise)\n\
402 -disccenter find rise of disc center (with option -rise)\n\
403 -discbottom find rise of disc bottom (with option -rise)\n\
404 -hindu hindu version of sunrise (with option -rise)\n\
405 specifications for heliacal events:\n\
406 -at[press,temp,rhum,visr]:\n\
407 pressure in hPa\n\
408 temperature in degrees Celsius\n\
409 relative humidity in %\n\
410 visual range, interpreted as follows:\n\
411 > 1 : meteorological range in km\n\
412 1>visr>0 : total atmospheric coefficient (ktot)\n\
413 = 0 : calculated from press, temp, rhum\n\
414 Default values are -at1013.25,15,40,0\n\
415 -obs[age,SN] age of observer and Snellen ratio\n\
416 Default values are -obs36,1\n\
417 -opt[age,SN,binocular,magn,diam,transm]\n\
418 age and SN as with -obs\n\
419 0 monocular or 1 binocular\n\
420 telescope magnification\n\
421 optical aperture in mm\n\
422 optical transmission\n\
423 Default values: -opt36,1,1,1,0,0 (naked eye)\n\
424 backward search:\n\
425 -bwd\n";
426 /* characters still available:
427 ijklruv
428 */
429 static char *infoplan = "\n\
430 Planet selection letters:\n\
431 planetary lists:\n\
432 d (default) main factors 0123456789mtABCcg\n\
433 p main factors as above, plus main asteroids DEFGHI\n\
434 h ficticious factors J..X\n\
435 a all factors\n\
436 (the letters above can only appear as a single letter)\n\n\
437 single body numbers/letters:\n\
438 0 Sun (character zero)\n\
439 1 Moon (character 1)\n\
440 2 Mercury\n\
441 3 Venus\n\
442 4 Mars\n\
443 5 Jupiter\n\
444 6 Saturn\n\
445 7 Uranus\n\
446 8 Neptune\n\
447 9 Pluto\n\
448 m mean lunar node\n\
449 t true lunar node\n\
450 n nutation\n\
451 o obliquity of ecliptic\n\
452 q delta t\n\
453 y time equation\n\
454 b ayanamsha\n\
455 A mean lunar apogee (Lilith, Black Moon) \n\
456 B osculating lunar apogee \n\
457 c intp. lunar apogee \n\
458 g intp. lunar perigee \n\
459 C Earth (in heliocentric or barycentric calculation)\n\
460 For planets Jupiter to Pluto the center of body (COB) can be\n\
461 calculated using the additional parameter -cob\n\
462 dwarf planets, plutoids\n\
463 F Ceres\n\
464 9 Pluto\n\
465 s -xs136199 Eris\n\
466 s -xs136472 Makemake\n\
467 s -xs136108 Haumea\n\
468 some minor planets:\n\
469 D Chiron\n\
470 E Pholus\n\
471 G Pallas \n\
472 H Juno \n\
473 I Vesta \n\
474 s minor planet, with MPC number given in -xs\n\
475 some planetary moons and center of body of a planet:\n\
476 v with moon number given in -xv:\n\
477 v -xv9501 Io/Jupiter:\n\
478 v -xv9599 Jupiter, center of body (COB):\n\
479 v -xv94.. Mars moons:\n\
480 v -xv95.. Jupiter moons and COB:\n\
481 v -xv96.. Saturn moons and COB:\n\
482 v -xv97.. Uranus moons and COB:\n\
483 v -xv98.. Neptune moons and COB:\n\
484 v -xv99.. Pluto moons and COB:\n\
485 The numbers of the moons are given here: \n\
486 https://www.astro.com/ftp/swisseph/ephe/sat/plmolist.txt\n\
487 fixed stars:\n\
488 f fixed star, with name or number given in -xf option\n\
489 f -xfSirius Sirius\n\
490 fictitious objects:\n\
491 J Cupido \n\
492 K Hades \n\
493 L Zeus \n\
494 M Kronos \n\
495 N Apollon \n\
496 O Admetos \n\
497 P Vulkanus \n\
498 Q Poseidon \n\
499 R Isis (Sevin) \n\
500 S Nibiru (Sitchin) \n\
501 T Harrington \n\
502 U Leverrier's Neptune\n\
503 V Adams' Neptune\n\
504 W Lowell's Pluto\n\
505 X Pickering's Pluto\n\
506 Y Vulcan\n\
507 Z White Moon\n\
508 w Waldemath's dark Moon\n\
509 z hypothetical body, with number given in -xz\n\
510 sidereal time:\n\
511 x sidereal time\n\
512 e print a line of labels\n\
513 \n";
514 /* characters still available
515 CcEeMmOoqWwz
516 */
517 static char *infoform = "\n\
518 Output format SEQ letters:\n\
519 In the standard setting five columns of coordinates are printed with\n\
520 the default format PLBRS. You can change the default by providing an\n\
521 option like -fCCCC where CCCC is your sequence of columns.\n\
522 The coding of the sequence is like this:\n\
523 y year\n\
524 Y year.fraction_of_year\n\
525 p planet index\n\
526 P planet name\n\
527 J absolute juldate\n\
528 T date formatted like 23.02.1992 \n\
529 t date formatted like 920223 for 1992 february 23\n\
530 L longitude in degree ddd mm'ss\"\n\
531 l longitude decimal\n\
532 Z longitude ddsignmm'ss\"\n\
533 S speed in longitude in degree ddd:mm:ss per day\n\
534 SS speed for all values specified in fmt\n\
535 s speed longitude decimal (degrees/day)\n\
536 ss speed for all values specified in fmt\n\
537 B latitude degree\n\
538 b latitude decimal\n\
539 R distance decimal in AU\n\
540 r distance decimal in AU, Moon in seconds parallax\n\
541 W distance decimal in light years\n\
542 w distance decimal in km\n\
543 q relative distance (1000=nearest, 0=furthest)\n\
544 A right ascension in hh:mm:ss\n\
545 a right ascension hours decimal\n\
546 m Meridian distance \n\
547 z Zenith distance \n\
548 D declination degree\n\
549 d declination decimal\n\
550 I azimuth degree\n\
551 i azimuth decimal\n\
552 H altitude degree\n\
553 h altitude decimal\n\
554 K altitude (with refraction) degree\n\
555 k altitude (with refraction) decimal\n\
556 G house position in degrees\n\
557 g house position in degrees decimal\n\
558 j house number 1.0 - 12.99999\n\
559 X x-, y-, and z-coordinates ecliptical\n\
560 x x-, y-, and z-coordinates equatorial\n\
561 U unit vector ecliptical\n\
562 u unit vector equatorial\n\
563 Q l, b, r, dl, db, dr, a, d, da, dd\n\
564 n nodes (mean): ascending/descending (Me - Ne); longitude decimal\n\
565 N nodes (osculating): ascending/descending, longitude; decimal\n\
566 f apsides (mean): perihelion, aphelion, second focal point; longitude dec.\n\
567 F apsides (osc.): perihelion, aphelion, second focal point; longitude dec.\n\
568 + phase angle\n\
569 - phase\n\
570 * elongation\n\
571 / apparent diameter of disc (without refraction)\n\
572 = magnitude\n";
573 static char *infoform2 = "\
574 v (reserved)\n\
575 V (reserved)\n\
576 ";
577 static char *infodate = "\n\
578 Date entry:\n\
579 In the interactive mode, when you are asked for a start date,\n\
580 you can enter data in one of the following formats:\n\
581 \n\
582 1.2.1991 three integers separated by a nondigit character for\n\
583 day month year. Dates are interpreted as Gregorian\n\
584 after 4.10.1582 and as Julian Calendar before.\n\
585 Time is always set to midnight (0 h).\n\
586 If the three letters jul are appended to the date,\n\
587 the Julian calendar is used even after 1582.\n\
588 If the four letters greg are appended to the date,\n\
589 the Gregorian calendar is used even before 1582.\n\
590 \n\
591 j2400123.67 the letter j followed by a real number, for\n\
592 the absolute Julian daynumber of the start date.\n\
593 Fraction .5 indicates midnight, fraction .0\n\
594 indicates noon, other times of the day can be\n\
595 chosen accordingly.\n\
596 \n\
597 <RETURN> repeat the last entry\n\
598 \n\
599 . stop the program\n\
600 \n\
601 +20 advance the date by 20 days\n\
602 \n\
603 -10 go back in time 10 days\n";
604 static char *infoexamp = "\n\
605 \n\
606 Examples:\n\
607 \n\
608 swetest -p2 -b1.12.1900 -n15 -s2\n\
609 ephemeris of Mercury (-p2) starting on 1 Dec 1900,\n\
610 15 positions (-n15) in two-day steps (-s2)\n\
611 \n\
612 swetest -p2 -b1.12.1900 -n15 -s2 -fTZ -roundsec -g, -head\n\
613 same, but output format = date and zodiacal position (-fTZ),\n\
614 separated by comma (-g,) and rounded to seconds (-roundsec),\n\
615 without header (-head).\n\
616 \n\
617 swetest -ps -xs433 -b1.12.1900\n\
618 position of asteroid 433 Eros (-ps -xs433)\n\
619 \n\
620 swetest -pf -xfAldebaran -b1.1.2000\n\
621 position of fixed star Aldebaran \n\
622 \n\
623 swetest -p1 -d0 -b1.12.1900 -n10 -fPTl -head\n\
624 angular distance of moon (-p1) from sun (-d0) for 10\n\
625 consecutive days (-n10).\n\
626 \n\
627 swetest -p6 -DD -b1.12.1900 -n100 -s5 -fPTZ -head -roundmin\n\
628 Midpoints between Saturn (-p6) and Chiron (-DD) for 100\n\
629 consecutive steps (-n100) with 5-day steps (-s5) with\n\
630 longitude in degree-sign format (-f..Z) rounded to minutes (-roundmin)\n\
631 \n\
632 swetest -b5.1.2002 -p -house12.05,49.50,K -ut12:30\n\
633 Koch houses for a location in Germany at a given date and time\n\
634 \n\
635 swetest -b1.1.2016 -g -fTlbR -p0123456789Dmte -hor -n366 -roundsec\n\
636 tabular ephemeris (all planets Sun - Pluto, Chiron, mean node, true node)\n\
637 in one horizontal row, tab-separated, for 366 days. For each planet\n\
638 list longitude, latitude and geocentric distance.\n";
639 /**************************************************************/
640
641 #include "swephexp.h" /* this includes "sweodef.h" */
642 #include "swephlib.h"
643 #include "sweph.h"
644
645 /*
646 * programmers warning: It looks much worse than it is!
647 * Originally swetest.c was a small and simple test program to test
648 * the main functions of the Swiss Ephemeris and to demonstrate
649 * its precision.
650 * It compiles on Unix, on MSDOS and as a non-GUI utility on 16-bit
651 * and 32-bit windows.
652 * This portability has forced us into some clumsy constructs, which
653 * end to hide the actual simplicity of the use of Swiss Ephemeris.
654 * For example, the mechanism implemented here in swetest.c to find
655 * the binary ephemeris files overrides the much simpler mechanism
656 * inside the SwissEph library. This was necessary because we wanted
657 * swetest.exe to run directly off the CDROM and search with some
658 * intelligence for ephemeris files already installed on a system.
659 */
660
661 #if MSDOS
662 # include <direct.h>
663 # include <dos.h>
664 # ifdef _MSC_VER
665 # include <sys\types.h>
666 # endif
667 #ifdef __MINGW32__
668 # include <sys/stat.h>
669 #else
670 # include <sys\stat.h>
671 #endif
672 # include <float.h>
673 #else
674 # ifdef MACOS
675 # include <console.h>
676 # else
677 # include <sys/stat.h>
678 # endif
679 #endif
680
681 #define J2000 2451545.0 /* 2000 January 1.5 */
682 #define square_sum(x) (x[0]*x[0]+x[1]*x[1]+x[2]*x[2])
683 #define SEFLG_EPHMASK (SEFLG_JPLEPH|SEFLG_SWIEPH|SEFLG_MOSEPH)
684
685 #define BIT_ROUND_SEC 1
686 #define BIT_ROUND_MIN 2
687 #define BIT_ZODIAC 4
688 #define BIT_LZEROES 8
689
690 #define BIT_TIME_LZEROES 8
691 #define BIT_TIME_LMT 16
692 #define BIT_TIME_LAT 32
693 #define BIT_ALLOW_361 64
694
695 #define PLSEL_D "0123456789mtA"
696 #define PLSEL_P "0123456789mtABCcgDEFGHI"
697 #define PLSEL_H "JKLMNOPQRSTUVWXYZw"
698 #define PLSEL_A "0123456789mtABCcgDEFGHIJKLMNOPQRSTUVWXYZw"
699
700 #define DIFF_DIFF 'd'
701 #define DIFF_GEOHEL 'h'
702 #define DIFF_MIDP 'D'
703 #define MODE_HOUSE 1
704 #define MODE_LABEL 2
705 #define MODE_AYANAMSA 4
706
707 #define SEARCH_RANGE_LUNAR_CYCLES 20000
708
709 #define LEN_SOUT 1000 // length of output string variable
710 #define SIND(x) sin((x) * DEGTORAD)
711 #define COSD(x) cos((x) * DEGTORAD)
712 #define ACOSD(x) (acos((x)) * RADTODEG)
713
714 static char se_pname[AS_MAXCH];
715 static char *zod_nam[] = {"ar", "ta", "ge", "cn", "le", "vi",
716 "li", "sc", "sa", "cp", "aq", "pi"};
717
718 static char star[AS_MAXCH] = "algol", star2[AS_MAXCH];
719 static char sastno[AS_MAXCH] = "433";
720 static char shyp[AS_MAXCH] = "1";
721 static char *dms(double x, int32 iflag);
722 static int make_ephemeris_path(char *argv0, char *ephepath);
723 static int letter_to_ipl(int letter);
724 static int print_line(int mode, AS_BOOL is_first, int sid_mode);
725 static int do_special_event(double tjd, int32 ipl, char *star, int32 special_event, int32 special_mode, double *geopos, double *datm, double *dobs, char *serr) ;
726 static int32 orbital_elements(double tjd_et, int32 ipl, int32 iflag, char *serr);
727 static char *hms_from_tjd(double x);
728 static void do_printf(char *info);
729 static char *hms(double x, int32 iflag);
730 static void remove_whitespace(char *s);
731 #if MSDOS
732 static int cut_str_any(char *s, char *cutlist, char *cpos[], int nmax);
733 #endif
734 static int32 call_swe_fixstar(char *star, double te, int32 iflag, double *x, char *serr);
735 static void jd_to_time_string(double jut, char *stimeout);
736 static char *our_strcpy(char *to, char *from);
737
738 /* globals shared between main() and print_line() */
739 static char *fmt = "PLBRS";
740 static char *gap = " ";
741 static double t, te, tut, jut = 0;
742 static int jmon, jday, jyear;
743 static int ipl = SE_SUN, ipldiff = SE_SUN, nhouses = 12;
744 static int iplctr = SE_SUN;
745 static char spnam[AS_MAXCH], spnam2[AS_MAXCH], serr[AS_MAXCH];
746 static char serr_save[AS_MAXCH], serr_warn[AS_MAXCH];
747 static int gregflag = SE_GREG_CAL;
748 static int diff_mode = 0;
749 static AS_BOOL use_dms = FALSE;
750 static AS_BOOL universal_time = FALSE;
751 static AS_BOOL universal_time_utc = FALSE;
752 static int32 round_flag = 0;
753 static int32 time_flag = 0;
754 static AS_BOOL short_output = FALSE;
755 static AS_BOOL list_hor = FALSE;
756 static int32 special_event = 0;
757 static int32 special_mode = 0;
758 static AS_BOOL do_orbital_elements = FALSE;
759 static AS_BOOL hel_using_AV = FALSE;
760 static AS_BOOL with_header = TRUE;
761 static AS_BOOL with_chart_link = FALSE;
762 static double x[6], x2[6], xequ[6], xcart[6], xcartq[6], xobl[6], xaz[6], xt[6], hpos, hpos2, hposj, armc, xsv[6];
763 static int hpos_meth = 0;
764 static double geopos[10];
765 static double attr[20], tret[20], datm[4], dobs[6];
766 static int32 iflag = 0, iflag2; /* external flag: helio, geo... */
767 static char *hs_nam[] =
768 {"undef", "Ascendant", "MC", "ARMC", "Vertex", "equat. Asc.","co-Asc. W.Koch", "co-Asc Munkasey", "Polar Asc."};
769 static int direction = 1;
770 static AS_BOOL direction_flag = FALSE;
771 static AS_BOOL step_in_minutes = FALSE;
772 static AS_BOOL step_in_seconds = FALSE;
773 static AS_BOOL step_in_years = FALSE;
774 static AS_BOOL step_in_months = FALSE;
775 static int32 helflag = 0;
776 static double tjd = 2415020.5;
777 static int32 nstep = 1, istep;
778 static int32 search_flag = 0;
779 static char sout[LEN_SOUT];
780 static int32 whicheph = SEFLG_SWIEPH;
781 static char *psp;
782 static int32 norefrac = 0;
783 static int32 disccenter = 0;
784 static int32 discbottom = 0;
785 static int32 hindu = 0;
786 /* for test of old models only */
787 static char *astro_models;
788 static int do_set_astro_models = FALSE;
789 static char smod[2000];
790 static AS_BOOL inut = FALSE; /* for Astrodienst internal feature */
791 static AS_BOOL have_gap_parameter = FALSE;
792 static AS_BOOL use_swe_fixstar2 = FALSE;
793 static AS_BOOL output_extra_prec = FALSE;
794 static AS_BOOL show_file_limit = FALSE;
795
796 #define SP_LUNAR_ECLIPSE 1
797 #define SP_SOLAR_ECLIPSE 2
798 #define SP_OCCULTATION 3
799 #define SP_RISE_SET 4
800 #define SP_MERIDIAN_TRANSIT 5
801 #define SP_HELIACAL 6
802
803 # define SP_MODE_HOW 2 /* an option for Lunar */
804 # define SP_MODE_LOCAL 8 /* an option for Solar */
805 # define SP_MODE_HOCAL 4096
806
807 # define ECL_LUN_PENUMBRAL 1 /* eclipse types for hocal list */
808 # define ECL_LUN_PARTIAL 2
809 # define ECL_LUN_TOTAL 3
810 # define ECL_SOL_PARTIAL 4
811 # define ECL_SOL_ANNULAR 5
812 # define ECL_SOL_TOTAL 6
813
main(int argc,char * argv[])814 int main(int argc, char *argv[])
815 {
816 char sdate_save[AS_MAXCH];
817 char s1[AS_MAXCH], s2[AS_MAXCH];
818 char *sp, *sp2;
819 char *spno;
820 char *plsel = PLSEL_D;
821 #if HPUNIX
822 char hostname[80];
823 #endif
824 int i, j, n, iflag_f = -1, iflgt;
825 int line_count, line_limit = 36525; // days in a century
826 double daya;
827 double top_long = 0.0; /* Greenwich UK */
828 double top_lat = 51.5;
829 double top_elev = 0;
830 AS_BOOL have_geopos = FALSE;
831 int ihsy = 'P';
832 AS_BOOL do_houses = FALSE;
833 char ephepath[AS_MAXCH];
834 char fname[AS_MAXCH];
835 char sdate[AS_MAXCH];
836 char *begindate = NULL;
837 char stimein[AS_MAXCH];
838 char stimeout[AS_MAXCH];
839 int32 iflgret;
840 AS_BOOL is_first = TRUE;
841 AS_BOOL with_glp = FALSE;
842 AS_BOOL with_header_always = FALSE;
843 AS_BOOL do_ayanamsa = FALSE;
844 AS_BOOL do_planeto_centric = FALSE;
845 double aya_t0 = 0, aya_val0 = 0;
846 AS_BOOL no_speed = FALSE;
847 int32 sid_mode = SE_SIDM_FAGAN_BRADLEY;
848 double t2, tstep = 1, thour = 0;
849 double delt;
850 double tid_acc = 0;
851 datm[0] = 1013.25; datm[1] = 15; datm[2] = 40; datm[3] = 0;
852 dobs[0] = 0; dobs[1] = 0;
853 dobs[2] = 0; dobs[3] = 0; dobs[4] = 0; dobs[5] = 0;
854 serr[0] = serr_save[0] = serr_warn[0] = sdate_save[0] = '\0';
855 # ifdef MACOS
856 argc = ccommand(&argv); /* display the arguments window */
857 # endif
858 *stimein = '\0';
859 strcpy(ephepath, "");
860 strcpy(fname, SE_FNAME_DFT);
861 for (i = 1; i < argc; i++) {
862 if (strncmp(argv[i], "-utc", 4) == 0) {
863 universal_time = TRUE;
864 universal_time_utc = TRUE;
865 if (strlen(argv[i]) > 4) {
866 strncpy(stimein, argv[i] + 4, 30);
867 stimein[30] = '\0';
868 }
869 } else if (strncmp(argv[i], "-ut", 3) == 0) {
870 universal_time = TRUE;
871 if (strlen(argv[i]) > 3) {
872 strncpy(stimein, argv[i] + 3, 30);
873 stimein[30] = '\0';
874 }
875 } else if (strncmp(argv[i], "-glp", 4) == 0) {
876 with_glp = TRUE;
877 } else if (strncmp(argv[i], "-hor", 4) == 0) {
878 list_hor = TRUE;
879 } else if (strncmp(argv[i], "-head", 5) == 0) {
880 with_header = FALSE;
881 } else if (strncmp(argv[i], "+head", 5) == 0) {
882 with_header_always = TRUE;
883 } else if (strcmp(argv[i], "-j2000") == 0) {
884 iflag |= SEFLG_J2000;
885 } else if (strcmp(argv[i], "-icrs") == 0) {
886 iflag |= SEFLG_ICRS;
887 } else if (strcmp(argv[i], "-cob") == 0) {
888 iflag |= SEFLG_CENTER_BODY;
889 } else if (strncmp(argv[i], "-ay", 3) == 0) {
890 do_ayanamsa = TRUE;
891 sid_mode = atol(argv[i]+3);
892 /*swe_set_sid_mode(sid_mode, 0, 0);*/
893 } else if (strncmp(argv[i], "-sidt0", 6) == 0) {
894 iflag |= SEFLG_SIDEREAL;
895 sid_mode = atol(argv[i]+6);
896 if (sid_mode == 0)
897 sid_mode = SE_SIDM_FAGAN_BRADLEY;
898 sid_mode |= SE_SIDBIT_ECL_T0;
899 /*swe_set_sid_mode(sid_mode, 0, 0);*/
900 } else if (strncmp(argv[i], "-sidsp", 6) == 0) {
901 iflag |= SEFLG_SIDEREAL;
902 sid_mode = atol(argv[i]+6);
903 if (sid_mode == 0)
904 sid_mode = SE_SIDM_FAGAN_BRADLEY;
905 sid_mode |= SE_SIDBIT_SSY_PLANE;
906 } else if (strncmp(argv[i], "-sidudef", 8) == 0) {
907 iflag |= SEFLG_SIDEREAL;
908 sid_mode = SE_SIDM_USER;
909 strcpy(s1, argv[i] + 8);
910 aya_t0 = atof(s1);
911 if ((sp = strchr(s1, ',')) != NULL) {
912 aya_val0 = atof(sp+1);
913 }
914 if (strstr(sp, "jdisut") != NULL) {
915 sid_mode |= SE_SIDBIT_USER_UT;
916 }
917 /*swe_set_sid_mode(sid_mode, 0, 0);*/
918 } else if (strncmp(argv[i], "-sidbit", 7) == 0) {
919 sid_mode |= atoi(argv[i]+7);
920 } else if (strncmp(argv[i], "-sid", 4) == 0) {
921 iflag |= SEFLG_SIDEREAL;
922 sid_mode = atol(argv[i]+4);
923 /*if (sid_mode > 0)
924 swe_set_sid_mode(sid_mode, 0, 0);*/
925 } else if (strcmp(argv[i], "-jplhora") == 0) {
926 iflag |= SEFLG_JPLHOR_APPROX;
927 } else if (strcmp(argv[i], "-tpm") == 0) {
928 iflag |= SEFLG_TEST_PLMOON;
929 } else if (strcmp(argv[i], "-jplhor") == 0) {
930 iflag |= SEFLG_JPLHOR;
931 } else if (strncmp(argv[i], "-j", 2) == 0) {
932 begindate = argv[i] + 1;
933 } else if (strncmp(argv[i], "-ejpl", 5) == 0) {
934 whicheph = SEFLG_JPLEPH;
935 if (*(argv[i]+5) != '\0') {
936 strncpy(fname, argv[i]+5, AS_MAXCH - 1);
937 fname[AS_MAXCH-1] = '\0';
938 }
939 } else if (strncmp(argv[i], "-edir", 5) == 0) {
940 if (*(argv[i]+5) != '\0') {
941 strncpy(ephepath, argv[i]+5, AS_MAXCH - 1);
942 ephepath[AS_MAXCH-1] = '\0';
943 }
944 } else if (strcmp(argv[i], "-eswe") == 0) {
945 whicheph = SEFLG_SWIEPH;
946 } else if (strcmp(argv[i], "-emos") == 0) {
947 whicheph = SEFLG_MOSEPH;
948 } else if (strncmp(argv[i], "-helflag", 8) == 0) {
949 helflag = atoi(argv[i]+8);
950 if (helflag >= SE_HELFLAG_AV)
951 hel_using_AV = TRUE;
952 } else if (strcmp(argv[i], "-hel") == 0) {
953 iflag |= SEFLG_HELCTR;
954 } else if (strcmp(argv[i], "-bary") == 0) {
955 iflag |= SEFLG_BARYCTR;
956 } else if (strncmp(argv[i], "-house", 6) == 0) {
957 sout[0] = '\0';
958 sout[1] = '\0';
959 sp = argv[i] + 6;
960 if (*sp == '[') sp++;
961 sscanf(sp, "%lf,%lf,%c", &top_long, &top_lat, sout);
962 top_elev = 0;
963 if (*sout) ihsy = sout[0];
964 do_houses = TRUE;
965 have_geopos = TRUE;
966 } else if (strncmp(argv[i], "-hsy", 4) == 0) {
967 ihsy = *(argv[i] + 4);
968 if (ihsy == '\0') ihsy = 'P';
969 if (strlen(argv[i]) > 5)
970 hpos_meth = atoi(argv[i] + 5);
971 have_geopos = TRUE;
972 } else if (strncmp(argv[i], "-topo", 5) == 0) {
973 iflag |= SEFLG_TOPOCTR;
974 sp = argv[i] + 5;
975 if (*sp == '[') sp++;
976 sscanf(sp, "%lf,%lf,%lf", &top_long, &top_lat, &top_elev);
977 have_geopos = TRUE;
978 } else if (strncmp(argv[i], "-geopos", 7) == 0) {
979 sp = argv[i] + 7;
980 if (*sp == '[') sp++;
981 sscanf(sp, "%lf,%lf,%lf", &top_long, &top_lat, &top_elev);
982 have_geopos = TRUE;
983 } else if (strcmp(argv[i], "-true") == 0) {
984 iflag |= SEFLG_TRUEPOS;
985 } else if (strcmp(argv[i], "-noaberr") == 0) {
986 iflag |= SEFLG_NOABERR;
987 } else if (strcmp(argv[i], "-nodefl") == 0) {
988 iflag |= SEFLG_NOGDEFL;
989 } else if (strcmp(argv[i], "-nonut") == 0) {
990 iflag |= SEFLG_NONUT;
991 } else if (strcmp(argv[i], "-speed3") == 0) {
992 iflag |= SEFLG_SPEED3;
993 } else if (strcmp(argv[i], "-speed") == 0) {
994 iflag |= SEFLG_SPEED;
995 } else if (strcmp(argv[i], "-nospeed") == 0) {
996 no_speed = TRUE;
997 } else if (strncmp(argv[i], "-testaa", 7) == 0) {
998 whicheph = SEFLG_JPLEPH;
999 strcpy(fname, SE_FNAME_DE200);
1000 if (strcmp(argv[i]+7, "95") == 0)
1001 begindate = "j2449975.5";
1002 if (strcmp(argv[i]+7, "96") == 0)
1003 begindate = "j2450442.5";
1004 if (strcmp(argv[i]+7, "97") == 0)
1005 begindate = "j2450482.5";
1006 fmt = "PADRu";
1007 universal_time = FALSE;
1008 plsel="3";
1009 } else if (strncmp(argv[i], "-lmt", 4) == 0) {
1010 universal_time = TRUE;
1011 time_flag |= BIT_TIME_LMT;
1012 if (strlen(argv[i]) > 4) {
1013 strncpy(stimein, argv[i] + 4, 30);
1014 stimein[30] = '\0';
1015 }
1016 } else if (strcmp(argv[i], "-lat") == 0) {
1017 universal_time = TRUE;
1018 time_flag |= BIT_TIME_LAT;
1019 } else if (strcmp(argv[i], "-lim") == 0) {
1020 show_file_limit = TRUE;
1021 } else if (strcmp(argv[i], "-clink") == 0) {
1022 with_chart_link = TRUE;
1023 } else if (strcmp(argv[i], "-lunecl") == 0) {
1024 special_event = SP_LUNAR_ECLIPSE;
1025 } else if (strcmp(argv[i], "-solecl") == 0) {
1026 special_event = SP_SOLAR_ECLIPSE;
1027 have_geopos = TRUE;
1028 } else if (strcmp(argv[i], "-short") == 0) {
1029 short_output = TRUE;
1030 } else if (strcmp(argv[i], "-occult") == 0) {
1031 special_event = SP_OCCULTATION;
1032 have_geopos = TRUE;
1033 } else if (strcmp(argv[i], "-ep") == 0) {
1034 output_extra_prec = TRUE;
1035 } else if (strcmp(argv[i], "-hocal") == 0) {
1036 /* used to create a listing for inclusion in hocal.c source code */
1037 special_mode |= SP_MODE_HOCAL;
1038 } else if (strcmp(argv[i], "-how") == 0) {
1039 special_mode |= SP_MODE_HOW;
1040 } else if (strcmp(argv[i], "-total") == 0) {
1041 search_flag |= SE_ECL_TOTAL;
1042 } else if (strcmp(argv[i], "-annular") == 0) {
1043 search_flag |= SE_ECL_ANNULAR;
1044 } else if (strcmp(argv[i], "-anntot") == 0) {
1045 search_flag |= SE_ECL_ANNULAR_TOTAL;
1046 } else if (strcmp(argv[i], "-partial") == 0) {
1047 search_flag |= SE_ECL_PARTIAL;
1048 } else if (strcmp(argv[i], "-penumbral") == 0) {
1049 search_flag |= SE_ECL_PENUMBRAL;
1050 } else if (strcmp(argv[i], "-noncentral") == 0) {
1051 search_flag &= ~SE_ECL_CENTRAL;
1052 search_flag |= SE_ECL_NONCENTRAL;
1053 } else if (strcmp(argv[i], "-central") == 0) {
1054 search_flag &= ~SE_ECL_NONCENTRAL;
1055 search_flag |= SE_ECL_CENTRAL;
1056 } else if (strcmp(argv[i], "-local") == 0) {
1057 special_mode |= SP_MODE_LOCAL;
1058 } else if (strcmp(argv[i], "-rise") == 0) {
1059 special_event = SP_RISE_SET;
1060 have_geopos = TRUE;
1061 } else if (strcmp(argv[i], "-norefrac") == 0) {
1062 norefrac = 1;
1063 } else if (strcmp(argv[i], "-disccenter") == 0) {
1064 disccenter = 1;
1065 } else if (strcmp(argv[i], "-hindu") == 0) {
1066 hindu = 1;
1067 norefrac = 1;
1068 disccenter = 1;
1069 } else if (strcmp(argv[i], "-discbottom") == 0) {
1070 discbottom = 1;
1071 } else if (strcmp(argv[i], "-metr") == 0) {
1072 special_event = SP_MERIDIAN_TRANSIT;
1073 have_geopos = TRUE;
1074 /* undocumented test feature */
1075 } else if (strncmp(argv[i], "-amod",5) == 0) {
1076 astro_models = argv[i] + 5;
1077 do_set_astro_models = TRUE;
1078 /* undocumented test feature */
1079 } else if (strncmp(argv[i], "-tidacc",7) == 0) {
1080 tid_acc = atof(argv[i] + 7);
1081 } else if (strncmp(argv[i], "-hev", 4) == 0) {
1082 special_event = SP_HELIACAL;
1083 search_flag = 0;
1084 sp = argv[i] + 4;
1085 if (*sp == '[') sp++;
1086 if (strlen(sp) > 0)
1087 search_flag = atoi(sp);
1088 have_geopos = TRUE;
1089 if (strstr(argv[i], "AV")) hel_using_AV = TRUE;
1090 } else if (strncmp(argv[i], "-at", 3) == 0) {
1091 sp = argv[i]+3;
1092 if (*sp == '[') sp++;
1093 j = 0;
1094 while (j < 4 && sp != NULL) {
1095 datm[j] = atof(sp);
1096 sp = strchr(sp, ',');
1097 if (sp != NULL) sp += 1;
1098 j++;
1099 }
1100 } else if (strncmp(argv[i], "-obs", 4) == 0) {
1101 sp = argv[i] + 4;
1102 if (*sp == '[') sp++;
1103 sscanf(sp, "%lf,%lf", &(dobs[0]), &(dobs[1]));
1104 } else if (strncmp(argv[i], "-opt", 4) == 0) {
1105 sp = argv[i] + 4;
1106 if (*sp == '[') sp++;
1107 sscanf(sp, "%lf,%lf,%lf,%lf,%lf,%lf", &(dobs[0]), &(dobs[1]), &(dobs[2]), &(dobs[3]), &(dobs[4]), &(dobs[5]));
1108 } else if (strcmp(argv[i], "-orbel") == 0) {
1109 do_orbital_elements = TRUE;
1110 } else if (strcmp(argv[i], "-bwd") == 0) {
1111 direction = -1;
1112 direction_flag = TRUE;
1113 } else if (strncmp(argv[i], "-pc", 3) == 0) {
1114 iplctr = atoi(argv[i]+3);
1115 do_planeto_centric = TRUE;
1116 } else if (strncmp(argv[i], "-p", 2) == 0) {
1117 spno = argv[i]+2;
1118 switch (*spno) {
1119 case 'd':
1120 /*
1121 case '\0':
1122 case ' ':
1123 */
1124 plsel = PLSEL_D; break;
1125 case 'p': plsel = PLSEL_P; break;
1126 case 'h': plsel = PLSEL_H; break;
1127 case 'a': plsel = PLSEL_A; break;
1128 default: plsel = spno;
1129 }
1130 } else if (strncmp(argv[i], "-xs", 3) == 0) {
1131 /* number of asteroid */
1132 strncpy(sastno, argv[i] + 3, AS_MAXCH - 1);
1133 sastno[AS_MAXCH-1] = '\0';
1134 } else if (strncmp(argv[i], "-xv", 3) == 0) {
1135 /* number of planetary moon */
1136 strncpy(sastno, argv[i] + 3, AS_MAXCH - 1);
1137 sastno[AS_MAXCH-1] = '\0';
1138 } else if (strncmp(argv[i], "-xf", 3) == 0) {
1139 /* name or number of fixed star */
1140 strncpy(star, argv[i] + 3, AS_MAXCH - 1);
1141 star[AS_MAXCH-1] = '\0';
1142 } else if (strncmp(argv[i], "-xz", 3) == 0) {
1143 /* number of hypothetical body */
1144 strncpy(shyp, argv[i] + 3, AS_MAXCH - 1);
1145 shyp[AS_MAXCH-1] = '\0';
1146 } else if (strncmp(argv[i], "-x", 2) == 0) {
1147 /* name or number of fixed star */
1148 strncpy(star, argv[i] + 2, AS_MAXCH - 1);
1149 star[AS_MAXCH-1] = '\0';
1150 } else if (strcmp(argv[i], "-nut") == 0) {
1151 inut = TRUE;
1152 } else if (strncmp(argv[i], "-n", 2) == 0) {
1153 nstep = atoi(argv[i]+2);
1154 if (nstep == 0)
1155 nstep = 20;
1156 } else if (strncmp(argv[i], "-i", 2) == 0) {
1157 iflag_f = atoi(argv[i]+2);
1158 if (iflag_f & SEFLG_XYZ)
1159 fmt = "PX";
1160 } else if (strcmp(argv[i], "-swefixstar2") == 0) {
1161 use_swe_fixstar2 = TRUE;
1162 } else if (strncmp(argv[i], "-s", 2) == 0) {
1163 tstep = atof(argv[i]+2);
1164 if (*(argv[i] + strlen(argv[i]) -1) == 'm')
1165 step_in_minutes = TRUE;
1166 if (*(argv[i] + strlen(argv[i]) -1) == 's')
1167 step_in_seconds = TRUE;
1168 if (*(argv[i] + strlen(argv[i]) -1) == 'y')
1169 step_in_years = TRUE;
1170 if (*(argv[i] + strlen(argv[i]) -1) == 'o') {
1171 step_in_minutes = FALSE;
1172 step_in_months = TRUE;
1173 }
1174 } else if (strncmp(argv[i], "-b", 2) == 0) {
1175 begindate = argv[i] + 2;
1176 } else if (strncmp(argv[i], "-f", 2) == 0) {
1177 fmt = argv[i] + 2;
1178 } else if (strncmp(argv[i], "-g", 2) == 0) {
1179 gap = argv[i] + 2;
1180 have_gap_parameter = TRUE;
1181 if (*gap == '\0') gap = "\t";
1182 } else if (strcmp(argv[i], "-dms") == 0) {
1183 use_dms = TRUE;
1184 } else if (strncmp(argv[i], "-d", 2) == 0
1185 || strncmp(argv[i], "-D", 2) == 0) {
1186 diff_mode = *(argv[i] + 1); /* 'd' or 'D' */
1187 sp = argv[i]+2;
1188 if (*(argv[i] + 2) == 'h') {
1189 sp++;
1190 diff_mode = 'h'; // diff helio to geo
1191 }
1192 ipldiff = letter_to_ipl((int) *sp);
1193 if (ipldiff <0) ipldiff = SE_SUN;
1194 swe_get_planet_name(ipldiff, spnam2);
1195 } else if (strcmp(argv[i], "-roundsec") == 0) {
1196 round_flag |= BIT_ROUND_SEC;
1197 } else if (strcmp(argv[i], "-roundmin") == 0) {
1198 round_flag |= BIT_ROUND_MIN;
1199 /*} else if (strncmp(argv[i], "-timeout", 8) == 0) {
1200 swe_set_timeout(atoi(argv[i]) + 8);*/
1201 } else if (strncmp(argv[i], "-t", 2) == 0) {
1202 if (strlen(argv[i]) > 2) {
1203 strncat(stimein, argv[i] + 2, 30);
1204 stimein[30] = '\0';
1205 }
1206 } else if (strncmp(argv[i], "-h", 2) == 0
1207 || strncmp(argv[i], "-?", 2) == 0) {
1208 sp = argv[i]+2;
1209 if (*sp == 'c' || *sp == '\0') {
1210 //char si0[strlen(infocmd0)+1]; // Microsoft Visual Studio does not like this
1211 char si0[2000];
1212 swe_version(sout);
1213 strcpy(si0, infocmd0);
1214 sp2 = strstr(si0, "Version:");
1215 if (sp2 != NULL && strlen(sp2) > 10 + strlen(sout))
1216 strcpy(sp2 + 9, sout);
1217 fputs(si0,stdout);
1218 fputs(infocmd1,stdout);
1219 fputs(infocmd2,stdout);
1220 fputs(infocmd3,stdout);
1221 fputs(infocmd4,stdout);
1222 fputs(infocmd5,stdout);
1223 fputs(infocmd6,stdout);
1224 }
1225 if (*sp == 'p' || *sp == '\0')
1226 fputs(infoplan,stdout);
1227 if (*sp == 'f' || *sp == '\0') {
1228 fputs(infoform,stdout);
1229 fputs(infoform2,stdout);
1230 }
1231 if (*sp == 'd' || *sp == '\0')
1232 fputs(infodate,stdout);
1233 if (*sp == 'e' || *sp == '\0')
1234 fputs(infoexamp,stdout);
1235 goto end_main;
1236 } else {
1237 strcpy(sout, "illegal option ");
1238 strncat(sout, argv[i], 100);
1239 sout[100] = '\0';
1240 strcat(sout, "\n");
1241 fputs(sout,stdout);
1242 exit(1);
1243 }
1244 }
1245 if (special_event == SP_OCCULTATION ||
1246 special_event == SP_RISE_SET ||
1247 special_event == SP_MERIDIAN_TRANSIT ||
1248 special_event == SP_HELIACAL
1249 ) {
1250 ipl = letter_to_ipl(*plsel);
1251 if (*plsel == 'f') {
1252 ipl = SE_FIXSTAR;
1253 } else {
1254 if (*plsel == 's')
1255 ipl = atoi(sastno) + SE_AST_OFFSET;
1256 *star = '\0';
1257 }
1258 if (special_event == SP_OCCULTATION && ipl == 1)
1259 ipl = 2; /* no occultation of moon by moon */
1260 }
1261 if (*stimein != '\0') {
1262 t = 0;
1263 if ((sp = strchr(stimein, ':')) != NULL) {
1264 if ((sp2 = strchr(sp + 1, ':')) != NULL) {
1265 t += atof(sp2 + 1) / 60.0;
1266 }
1267 t += atoi(sp + 1);
1268 t /= 60.0;
1269 }
1270 if (atoi(stimein) < 0)
1271 t = -t;
1272 t += atoi(stimein);
1273 //t += 0.0000000001;
1274 thour = t;
1275 }
1276 #if HPUNIX
1277 gethostname (hostname, 80);
1278 if (strstr(hostname, "as10") != NULL)
1279 line_limit = 1000;
1280 #endif
1281 #if MSDOS
1282 SetConsoleOutputCP(65001); // set console to utf-8,
1283 // works only from Windows Vista upwards, not on XP.
1284 #endif
1285 if (with_header) {
1286 for (i = 0; i < argc; i++) {
1287 fputs(argv[i],stdout);
1288 printf(" ");
1289 }
1290 }
1291 iflag = (iflag & ~SEFLG_EPHMASK) | whicheph;
1292 if (strpbrk(fmt, "SsQ") != NULL && !(iflag & SEFLG_SPEED3) && !no_speed)
1293 iflag |= SEFLG_SPEED;
1294 if (*ephepath == '\0') {
1295 if (make_ephemeris_path(argv[0], ephepath) == ERR) {
1296 iflag = (iflag & ~SEFLG_EPHMASK) | SEFLG_MOSEPH;
1297 whicheph = SEFLG_MOSEPH;
1298 }
1299 }
1300 if (whicheph != SEFLG_MOSEPH)
1301 swe_set_ephe_path(ephepath);
1302 if (whicheph & SEFLG_JPLEPH)
1303 swe_set_jpl_file(fname);
1304 /* the following is only a test feature */
1305 if (do_set_astro_models) {
1306 swe_set_astro_models(astro_models, iflag); /* secret test feature for dieter */
1307 swe_get_astro_models(astro_models, smod, iflag);
1308 }
1309 #if 1
1310 if (inut) /* Astrodienst internal feature */
1311 swe_set_interpolate_nut(TRUE);
1312 #endif
1313 if ((iflag & SEFLG_SIDEREAL) || do_ayanamsa) {
1314 if (sid_mode & SE_SIDM_USER)
1315 swe_set_sid_mode(sid_mode, aya_t0, aya_val0);
1316 else
1317 swe_set_sid_mode(sid_mode, 0, 0);
1318 }
1319 geopos[0] = top_long;
1320 geopos[1] = top_lat;
1321 geopos[2] = top_elev;
1322 swe_set_topo(top_long, top_lat, top_elev);
1323 if (tid_acc != 0)
1324 swe_set_tid_acc(tid_acc);
1325 serr[0] = serr_save[0] = serr_warn[0] = '\0';
1326 while (TRUE) {
1327 if (begindate == NULL) {
1328 printf("\nDate ?");
1329 sdate[0] = '\0';
1330 if( !fgets(sdate, AS_MAXCH, stdin) ) goto end_main;
1331 } else {
1332 strncpy(sdate, begindate, AS_MAXCH-1);
1333 sdate[AS_MAXCH-1] = '\0';
1334 begindate = "."; /* to exit afterwards */
1335 }
1336 if (strcmp(sdate, "-bary") == 0) {
1337 iflag = iflag & ~SEFLG_HELCTR;
1338 iflag |= SEFLG_BARYCTR;
1339 *sdate = '\0';
1340 } else if (strcmp(sdate, "-hel") == 0) {
1341 iflag = iflag & ~SEFLG_BARYCTR;
1342 iflag |= SEFLG_HELCTR;
1343 *sdate = '\0';
1344 } else if (strcmp(sdate, "-geo") == 0) {
1345 iflag = iflag & ~SEFLG_BARYCTR;
1346 iflag = iflag & ~SEFLG_HELCTR;
1347 *sdate = '\0';
1348 } else if (strcmp(sdate, "-ejpl") == 0) {
1349 iflag &= ~SEFLG_EPHMASK;
1350 iflag |= SEFLG_JPLEPH;
1351 *sdate = '\0';
1352 } else if (strcmp(sdate, "-eswe") == 0) {
1353 iflag &= ~SEFLG_EPHMASK;
1354 iflag |= SEFLG_SWIEPH;
1355 *sdate = '\0';
1356 } else if (strcmp(sdate, "-emos") == 0) {
1357 iflag &= ~SEFLG_EPHMASK;
1358 iflag |= SEFLG_MOSEPH;
1359 *sdate = '\0';
1360 } else if (strncmp(sdate, "-xs",3) == 0) {
1361 /* number of asteroid */
1362 strcpy(sastno, sdate + 3);
1363 *sdate = '\0';
1364 }
1365 sp = sdate;
1366 if (*sp == '.') {
1367 goto end_main;
1368 } else if (*sp == '\0' || *sp == '\n' || *sp == '\r') {
1369 strcpy(sdate, sdate_save);
1370 } else {
1371 strcpy(sdate_save, sdate);
1372 }
1373 if (*sdate == '\0') {
1374 sprintf(sdate, "j%f", tjd);
1375 }
1376 if (*sp == 'j') { /* it's a day number */
1377 if ((sp2 = strchr(sp, ',')) != NULL)
1378 *sp2 = '.';
1379 sscanf(sp+1,"%lf", &tjd);
1380 if (tjd < 2299160.5)
1381 gregflag = SE_JUL_CAL;
1382 else
1383 gregflag = SE_GREG_CAL;
1384 if (strstr(sp, "jul") != NULL)
1385 gregflag = SE_JUL_CAL;
1386 else if (strstr(sp, "greg") != NULL)
1387 gregflag = SE_GREG_CAL;
1388 swe_revjul(tjd, gregflag, &jyear, &jmon, &jday, &jut);
1389 } else if (*sp == '+') {
1390 n = atoi(sp);
1391 if (n == 0) n = 1;
1392 tjd += n;
1393 swe_revjul(tjd, gregflag, &jyear, &jmon, &jday, &jut);
1394 } else if (*sp == '-') {
1395 n = atoi(sp);
1396 if (n == 0) n = -1;
1397 tjd += n;
1398 swe_revjul(tjd, gregflag, &jyear, &jmon, &jday, &jut);
1399 } else {
1400 if (sscanf (sp, "%d%*c%d%*c%d", &jday,&jmon,&jyear) < 1) exit(1);
1401 if ((int32) jyear * 10000L + (int32) jmon * 100L + (int32) jday < 15821015L)
1402 gregflag = SE_JUL_CAL;
1403 else
1404 gregflag = SE_GREG_CAL;
1405 if (strstr(sp, "jul") != NULL)
1406 gregflag = SE_JUL_CAL;
1407 else if (strstr(sp, "greg") != NULL)
1408 gregflag = SE_GREG_CAL;
1409 jut = 0;
1410 if (universal_time_utc) {
1411 int ih = 0, im = 0;
1412 double ds = 0.0;
1413 if (*stimein != '\0') {
1414 sscanf(stimein, "%d:%d:%lf", &ih, &im, &ds);
1415 }
1416 if (swe_utc_to_jd(jyear,jmon,jday, ih, im, ds, gregflag, tret, serr) == ERR) {
1417 printf(" error in swe_utc_to_jd(): %s\n", serr);
1418 exit(-1);
1419 }
1420 tjd = tret[1];
1421 } else {
1422 tjd = swe_julday(jyear,jmon,jday,jut,gregflag);
1423 tjd += thour / 24.0;
1424 }
1425 }
1426 if (special_event > 0) {
1427 do_special_event(tjd, ipl, star, special_event, special_mode, geopos, datm, dobs, serr) ;
1428 swe_close();
1429 return OK;
1430 }
1431 line_count = 0;
1432 for (t = tjd, istep = 1; istep <= nstep; t += tstep, istep++) {
1433 if (step_in_minutes)
1434 t = tjd + (istep -1) * tstep / 1440;
1435 if (step_in_seconds)
1436 t = tjd + (istep -1) * tstep / 86400;
1437 if (step_in_years) {
1438 swe_revjul(tjd, gregflag, &jyear, &jmon, &jday, &jut);
1439 t = swe_julday(jyear + (istep - 1) * (int) tstep, jmon, jday, jut, gregflag);
1440 }
1441 if (step_in_months) {
1442 swe_revjul(tjd, gregflag, &jyear, &jmon, &jday, &jut);
1443 jmon += (istep - 1) * (int) tstep;
1444 jyear += (int) ((jmon - 1) / 12);
1445 jmon = ((jmon - 1) % 12) + 1;
1446 t = swe_julday(jyear, jmon, jday, jut, gregflag);
1447 }
1448 if (t < 2299160.5)
1449 gregflag = SE_JUL_CAL;
1450 else
1451 gregflag = SE_GREG_CAL;
1452 if (strstr(sdate, "jul") != NULL)
1453 gregflag = SE_JUL_CAL;
1454 else if (strstr(sdate, "greg") != NULL)
1455 gregflag = SE_GREG_CAL;
1456 delt = swe_deltat_ex(t, iflag, serr);
1457 if (!universal_time) {
1458 delt = swe_deltat_ex(t - delt, iflag, serr);
1459 }
1460 t2 = t;
1461 // output line:
1462 // "date (dmy) 4.6.2017 greg. 2:07:00 TT version 2.07.02"
1463 swe_revjul(t2, gregflag, &jyear, &jmon, &jday, &jut);
1464 if (with_header) {
1465 #ifndef NO_SWE_GLP // -DNO_SWE_GLP to suppress this function
1466 if (with_glp) {
1467 swe_get_library_path(sout);
1468 printf("\npath: %s", sout);
1469 }
1470 #endif
1471 printf("\ndate (dmy) %d.%d.%d", jday, jmon, jyear);
1472 if (gregflag)
1473 printf(" greg.");
1474 else
1475 printf(" jul.");
1476 jd_to_time_string(jut, stimeout);
1477 printf(stimeout);
1478 if (universal_time) {
1479 if (time_flag & BIT_TIME_LMT)
1480 printf(" LMT");
1481 else
1482 printf(" UT");
1483 } else {
1484 printf(" TT");
1485 }
1486 printf("\t\tversion %s", swe_version(sout));
1487 }
1488 if (universal_time) {
1489 // "LMT: 2457908.588194444"
1490 if (time_flag & BIT_TIME_LMT) {
1491 if (with_header) {
1492 printf("\nLMT: %.9f", t);
1493 t -= geopos[0] / 15.0 / 24.0;
1494 }
1495 }
1496 // "UT: 2457908.565972222 delta t: 68.761612 sec"
1497 if (with_header) {
1498 printf("\nUT: %.9f", t);
1499 printf(" delta t: %f sec", delt * 86400.0);
1500 }
1501 te = t + delt;
1502 tut = t;
1503 } else {
1504 te = t;
1505 tut = t - delt;
1506 // "UT: 2457908.565972222 delta t: 68.761612 sec"
1507 if (with_header) {
1508 printf("\nUT: %.9f", tut);
1509 printf(" delta t: %f sec", delt * 86400.0);
1510 }
1511 }
1512 iflgret = swe_calc(te, SE_ECL_NUT, iflag, xobl, serr);
1513 if (with_header) {
1514 // "TT: 2457908.566768074
1515 printf("\nTT: %.9f", te);
1516 // "ayanamsa = 24° 5'51.6509 (Lahiri)"
1517 if (iflag & SEFLG_SIDEREAL) {
1518 if (swe_get_ayanamsa_ex(te, iflag, &daya, serr) == ERR) {
1519 printf(" error in swe_get_ayanamsa_ex(): %s\n", serr);
1520 exit(1);
1521 }
1522 printf(" ayanamsa = %s (%s)", dms(daya, round_flag), swe_get_ayanamsa_name(sid_mode));
1523 }
1524 // "geo. long 8.000000, lat 47.000000, alt 0.000000"
1525 if (have_geopos) {
1526 printf("\ngeo. long %f, lat %f, alt %f", geopos[0], geopos[1], geopos[2]);
1527 }
1528 if (iflag_f >=0)
1529 iflag = iflag_f;
1530 if (strchr(plsel, 'o') == NULL) {
1531 if (iflag & (SEFLG_NONUT | SEFLG_SIDEREAL)) {
1532 printf("\n%-15s %s", "Epsilon (m)", dms(xobl[0],round_flag));
1533 } else {
1534 printf("\n%-15s %s%s", "Epsilon (t/m)", dms(xobl[0],round_flag), gap);
1535 printf("%s", dms(xobl[1],round_flag));
1536 }
1537 }
1538 if (strchr(plsel, 'n') == NULL && !(iflag & (SEFLG_NONUT | SEFLG_SIDEREAL))) {
1539 fputs("\nNutation ", stdout);
1540 fputs(dms(xobl[2], round_flag), stdout);
1541 fputs(gap, stdout);
1542 fputs(dms(xobl[3], round_flag), stdout);
1543 }
1544 printf("\n");
1545 if (do_houses) {
1546 char *shsy = swe_house_name(ihsy);
1547 if (!universal_time) {
1548 do_houses = FALSE;
1549 printf("option -house requires option -ut for Universal Time\n");
1550 } else {
1551 strcpy(s1, dms(top_long, round_flag));
1552 strcpy(s2, dms(top_lat, round_flag));
1553 printf("Houses system %c (%s) for long=%s, lat=%s\n", ihsy, shsy, s1, s2);
1554 }
1555 }
1556 }
1557 if (with_header && !with_header_always)
1558 with_header = FALSE;
1559 if (do_ayanamsa) {
1560 if (swe_get_ayanamsa_ex(te, iflag, &daya, serr) == ERR) {
1561 printf(" error in swe_get_ayanamsa_ex(): %s\n", serr);
1562 exit(1);
1563 }
1564 x[0] = daya;
1565 print_line(MODE_AYANAMSA, TRUE, sid_mode);
1566 continue;
1567 }
1568 if (t == tjd && strchr(plsel, 'e')) {
1569 if (list_hor) {
1570 is_first = TRUE;
1571 for (psp = plsel; *psp != '\0'; psp++) {
1572 if (*psp == 'e') continue;
1573 ipl = letter_to_ipl((int) *psp);
1574 *spnam = '\0';
1575 if (ipl >= SE_SUN && ipl <= SE_VESTA)
1576 swe_get_planet_name(ipl, spnam);
1577 print_line(MODE_LABEL, is_first, 0);
1578 is_first = FALSE;
1579 }
1580 printf("\n");
1581 } else {
1582 print_line(MODE_LABEL, TRUE, 0);
1583 }
1584 }
1585 is_first = TRUE;
1586 for (psp = plsel; *psp != '\0'; psp++) {
1587 if (*psp == 'e') continue;
1588 ipl = letter_to_ipl((int) *psp);
1589 if (ipl == -2) {
1590 printf("illegal parameter -p%s\n", plsel);
1591 exit(1);
1592 }
1593 if (*psp == 'f') // fixed star
1594 ipl = SE_FIXSTAR;
1595 else if (*psp == 's') // asteroid
1596 ipl = atoi(sastno) + 10000;
1597 else if (*psp == 'v') // planetary moon
1598 ipl = atoi(sastno);
1599 else if (*psp == 'z') // fictitious object
1600 ipl = atoi(shyp) + SE_FICT_OFFSET_1;
1601 if (iflag & SEFLG_HELCTR) {
1602 if (ipl == SE_SUN
1603 || ipl == SE_MEAN_NODE || ipl == SE_TRUE_NODE
1604 || ipl == SE_MEAN_APOG || ipl == SE_OSCU_APOG)
1605 continue;
1606 } else if (iflag & SEFLG_BARYCTR) {
1607 if (ipl == SE_MEAN_NODE || ipl == SE_TRUE_NODE
1608 || ipl == SE_MEAN_APOG || ipl == SE_OSCU_APOG)
1609 continue;
1610 } else { /* geocentric */
1611 if (ipl == SE_EARTH && !do_orbital_elements)
1612 continue;
1613 }
1614 /* ecliptic position */
1615 if (iflag_f >=0)
1616 iflag = iflag_f;
1617 if (ipl == SE_FIXSTAR) {
1618 iflgret = call_swe_fixstar(star, te, iflag, x, serr);
1619 /* magnitude, etc. */
1620 if (iflgret != ERR && strpbrk(fmt, "=") != NULL) {
1621 double mag;
1622 iflgret = swe_fixstar_mag(star, &mag, serr);
1623 attr[4] = mag;
1624 }
1625 strcpy(se_pname, star);
1626 } else if (do_planeto_centric) {
1627 iflgret = swe_calc_pctr(te, ipl, iplctr, iflag, x, serr);
1628 swe_get_planet_name(ipl, se_pname);
1629 } else {
1630 iflgret = swe_calc(te, ipl, iflag, x, serr);
1631 /* phase, magnitude, etc. */
1632 if (iflgret != ERR && strpbrk(fmt, "+-*/=") != NULL)
1633 iflgret = swe_pheno(te, ipl, iflag, attr, serr);
1634 swe_get_planet_name(ipl, se_pname);
1635 if (show_file_limit && ipl > SE_AST_OFFSET) {
1636 const char *fnam;
1637 char sbeg[40], send[40];
1638 double tfstart, tfend;
1639 int denum;
1640 fnam = swe_get_current_file_data(3, &tfstart, &tfend, &denum);
1641 if (fnam != NULL) {
1642 swe_revjul(tfstart, gregflag, &jyear, &jmon, &jday, &jut);
1643 sprintf(sbeg, "%d.%02d.%04d", jday, jmon, jyear);
1644 swe_revjul(tfend, gregflag, &jyear, &jmon, &jday, &jut);
1645 sprintf(send, "%d.%02d.%04d", jday, jmon, jyear);
1646 printf("range %s: %.1lf = %s to %.1lf = %s de=%d\n", fnam, tfstart, sbeg, tfend, send, denum);
1647 show_file_limit = FALSE;
1648 }
1649 }
1650 }
1651 if (*psp == 'q') {/* delta t */
1652 x[0] = swe_deltat_ex(tut, iflag, serr) * 86400;
1653 x[1] = x[2] = x[3] = 0;
1654 x[1] = x[0] / 3600.0; // to hours
1655 strcpy(se_pname, "Delta T");
1656 }
1657 if (*psp == 'x') {/* sidereal time */
1658 x[0] = swe_degnorm(swe_sidtime(tut) * 15 + geopos[0]);
1659 x[1] = x[2] = x[3] = 0;
1660 strcpy(se_pname, "Sidereal Time");
1661 }
1662 if (*psp == 'o') {/* ecliptic is wanted, remove nutation */
1663 x[2] = x[3] = 0;
1664 strcpy(se_pname, "Ecl. Obl.");
1665 }
1666 if (*psp == 'n') {/* nutation is wanted, remove ecliptic */
1667 x[0] = x[2];
1668 x[1] = x[3];
1669 x[2] = x[3] = 0;
1670 strcpy(se_pname, "Nutation");
1671 }
1672 if (*psp == 'y') {/* time equation */
1673 iflgret = swe_time_equ(tut, &(x[0]), serr);
1674 x[0] *= 86400; /* in seconds */;
1675 x[1] = x[2] = x[3] = 0;
1676 strcpy(se_pname, "Time Equ.");
1677 }
1678 if (*psp == 'b') {/* ayanamsha */
1679 if (swe_get_ayanamsa_ex(te, iflag, &(x[0]), serr) == ERR) {
1680 printf(" error in swe_get_ayanamsa_ex(): %s\n", serr);
1681 iflgret = -1;
1682 }
1683 x[1] = 0;
1684 strcpy(se_pname, "Ayanamsha");
1685 }
1686 if (iflgret < 0) {
1687 if (strcmp(serr, serr_save) != 0
1688 && (ipl == SE_SUN || ipl == SE_MOON || ipl <= SE_PLUTO
1689 || ipl == SE_MEAN_NODE || ipl == SE_TRUE_NODE
1690 || ipl == SE_CERES || ipl == SE_PALLAS || ipl == SE_JUNO || ipl == SE_VESTA
1691 || ipl == SE_CHIRON || ipl == SE_PHOLUS || ipl == SE_CUPIDO
1692 || ipl >= SE_PLMOON_OFFSET
1693 || ipl >= SE_AST_OFFSET || ipl == SE_FIXSTAR
1694 || *psp == 'y')) {
1695 fputs("error: ", stdout);
1696 fputs(serr, stdout);
1697 fputs("\n", stdout);
1698 }
1699 strcpy(serr_save, serr);
1700 } else if (*serr != '\0' && *serr_warn == '\0') {
1701 if (strstr(serr, "'seorbel.txt' not found") == NULL)
1702 strcpy(serr_warn, serr);
1703 }
1704 if (diff_mode) {
1705 iflgret = swe_calc(te, ipldiff, iflag, x2, serr);
1706 if (diff_mode == DIFF_GEOHEL)
1707 iflgret = swe_calc(te, ipldiff, iflag|SEFLG_HELCTR, x2, serr);
1708 if (iflgret < 0) {
1709 fputs("error: ", stdout);
1710 fputs(serr, stdout);
1711 fputs("\n", stdout);
1712 }
1713 if (diff_mode == DIFF_DIFF || diff_mode == DIFF_GEOHEL) {
1714 for (i = 1; i < 6; i++)
1715 x[i] -= x2[i];
1716 if ((iflag & SEFLG_RADIANS) == 0)
1717 x[0] = swe_difdeg2n(x[0], x2[0]);
1718 else
1719 x[0] = swe_difrad2n(x[0], x2[0]);
1720 } else { /* DIFF_MIDP */
1721 for (i = 1; i < 6; i++)
1722 x[i] = (x[i] + x2[i]) / 2;
1723 if ((iflag & SEFLG_RADIANS) == 0)
1724 x[0] = swe_deg_midp(x[0], x2[0]);
1725 else
1726 x[0] = swe_rad_midp(x[0], x2[0]);
1727 }
1728 }
1729 /* equator position */
1730 if (strpbrk(fmt, "aADdQmz") != NULL) {
1731 iflag2 = iflag | SEFLG_EQUATORIAL;
1732 if (ipl == SE_FIXSTAR) {
1733 iflgret = call_swe_fixstar(star, te, iflag2, xequ, serr);
1734 } else if (do_planeto_centric) {
1735 iflgret = swe_calc_pctr(te, ipl, iplctr, iflag2, xequ, serr);
1736 } else {
1737 iflgret = swe_calc(te, ipl, iflag2, xequ, serr);
1738 }
1739 if (diff_mode) {
1740 iflgret = swe_calc(te, ipldiff, iflag2, x2, serr);
1741 if (diff_mode == DIFF_DIFF || diff_mode == DIFF_GEOHEL) {
1742 if (diff_mode == DIFF_GEOHEL)
1743 iflgret = swe_calc(te, ipldiff, iflag2|SEFLG_HELCTR, x2, serr);
1744 for (i = 1; i < 6; i++)
1745 xequ[i] -= x2[i];
1746 if ((iflag & SEFLG_RADIANS) == 0)
1747 xequ[0] = swe_difdeg2n(xequ[0], x2[0]);
1748 else
1749 xequ[0] = swe_difrad2n(xequ[0], x2[0]);
1750 } else { /* DIFF_MIDP */
1751 for (i = 1; i < 6; i++)
1752 xequ[i] = (xequ[i] + x2[i]) / 2;
1753 if ((iflag & SEFLG_RADIANS) == 0)
1754 xequ[0] = swe_deg_midp(xequ[0], x2[0]);
1755 else
1756 xequ[0] = swe_rad_midp(xequ[0], x2[0]);
1757 }
1758 }
1759 }
1760 /* azimuth and height */
1761 if (strpbrk(fmt, "IiHhKk") != NULL) {
1762 /* first, get topocentric equatorial positions */
1763 iflgt = whicheph | SEFLG_EQUATORIAL | SEFLG_TOPOCTR;
1764 if (ipl == SE_FIXSTAR)
1765 iflgret = call_swe_fixstar(star, te, iflgt, xt, serr);
1766 else
1767 iflgret = swe_calc(te, ipl, iflgt, xt, serr);
1768 /* to azimuth/height */
1769 /* atmospheric pressure "0" has the effect that a value
1770 * of 1013.25 mbar is assumed at 0 m above sea level.
1771 * If the altitude of the observer is given (in geopos[2])
1772 * pressure is estimated according to that */
1773 swe_azalt(tut, SE_EQU2HOR, geopos, datm[0], datm[1], xt, xaz);
1774 if (diff_mode) {
1775 iflgret = swe_calc(te, ipldiff, iflgt, xt, serr);
1776 swe_azalt(tut, SE_EQU2HOR, geopos, datm[0], datm[1], xt, x2);
1777 if (diff_mode == DIFF_DIFF || diff_mode == DIFF_GEOHEL) {
1778 if (diff_mode == DIFF_GEOHEL) { // makes little sense for a heliocentric
1779 iflgret = swe_calc(te, ipldiff, iflgt|SEFLG_HELCTR, xt, serr);
1780 swe_azalt(tut, SE_EQU2HOR, geopos, datm[0], datm[1], xt, x2);
1781 }
1782 for (i = 1; i < 3; i++)
1783 xaz[i] -= x2[i];
1784 if ((iflag & SEFLG_RADIANS) == 0)
1785 xaz[0] = swe_difdeg2n(xaz[0], x2[0]);
1786 else
1787 xaz[0] = swe_difrad2n(xaz[0], x2[0]);
1788 } else { /* DIFF_MIDP */
1789 for (i = 1; i < 3; i++)
1790 xaz[i] = (xaz[i] + x2[i]) / 2;
1791 if ((iflag & SEFLG_RADIANS) == 0)
1792 xaz[0] = swe_deg_midp(xaz[0], x2[0]);
1793 else
1794 xaz[0] = swe_rad_midp(xaz[0], x2[0]);
1795 }
1796 }
1797 }
1798 /* ecliptic cartesian position */
1799 if (strpbrk(fmt, "XU") != NULL) {
1800 iflag2 = iflag | SEFLG_XYZ;
1801 if (ipl == SE_FIXSTAR)
1802 iflgret = call_swe_fixstar(star, te, iflag2, xcart, serr);
1803 else
1804 iflgret = swe_calc(te, ipl, iflag2, xcart, serr);
1805 if (diff_mode) {
1806 iflgret = swe_calc(te, ipldiff, iflag2, x2, serr);
1807 if (diff_mode == DIFF_DIFF || diff_mode == DIFF_GEOHEL) {
1808 if (diff_mode == DIFF_GEOHEL)
1809 iflgret = swe_calc(te, ipldiff, iflag2|SEFLG_HELCTR, x2, serr);
1810 for (i = 0; i < 6; i++)
1811 xcart[i] -= x2[i];
1812 } else {
1813 xcart[i] = (xcart[i] + x2[i]) / 2;
1814 }
1815 }
1816 }
1817 /* equator cartesian position */
1818 if (strpbrk(fmt, "xu") != NULL) {
1819 iflag2 = iflag | SEFLG_XYZ | SEFLG_EQUATORIAL;
1820 if (ipl == SE_FIXSTAR)
1821 iflgret = call_swe_fixstar(star, te, iflag2, xcartq, serr);
1822 else
1823 iflgret = swe_calc(te, ipl, iflag2, xcartq, serr);
1824 if (diff_mode) {
1825 iflgret = swe_calc(te, ipldiff, iflag2, x2, serr);
1826 if (diff_mode == DIFF_DIFF || diff_mode == DIFF_GEOHEL) {
1827 if (diff_mode == DIFF_GEOHEL)
1828 iflgret = swe_calc(te, ipldiff, iflag2|SEFLG_HELCTR, x2, serr);
1829 for (i = 0; i < 6; i++)
1830 xcartq[i] -= x2[i];
1831 } else {
1832 xcartq[i] = (xcart[i] + x2[i]) / 2;
1833 }
1834 }
1835 }
1836 /* house position */
1837 if (strpbrk(fmt, "gGjzm") != NULL) {
1838 armc = swe_degnorm(swe_sidtime(tut) * 15 + geopos[0]);
1839 for (i = 0; i < 6; i++)
1840 xsv[i] = x[i];
1841 if (hpos_meth == 1)
1842 xsv[1] = 0;
1843 if (ipl == SE_FIXSTAR)
1844 strcpy(star2, star);
1845 else
1846 *star2 = '\0';
1847 if (hpos_meth >= 2 && toupper(ihsy) == 'G') {
1848 swe_gauquelin_sector(tut, ipl, star2, iflag, hpos_meth, geopos, 0, 0, &hposj, serr);
1849 } else {
1850 double cusp[100];
1851 if (ihsy == 'i' || ihsy == 'I') {
1852 iflgret = swe_houses_ex(t,iflag, top_lat, top_long, ihsy, cusp, cusp+13);
1853 }
1854 hposj = swe_house_pos(armc, geopos[1], xobl[0], ihsy, xsv, serr);
1855 }
1856 if (toupper(ihsy) == 'G')
1857 hpos = (hposj - 1) * 10;
1858 else
1859 hpos = (hposj - 1) * 30;
1860 if (diff_mode) {
1861 for (i = 0; i < 6; i++)
1862 xsv[i] = x2[i];
1863 if (hpos_meth == 1)
1864 xsv[1] = 0;
1865 hpos2 = swe_house_pos(armc, geopos[1], xobl[0], ihsy, xsv, serr);
1866 if (toupper(ihsy) == 'G')
1867 hpos2 = (hpos2 - 1) * 10;
1868 else
1869 hpos2 = (hpos2 - 1) * 30;
1870 if (diff_mode == DIFF_DIFF || diff_mode == DIFF_GEOHEL) {
1871 if ((iflag & SEFLG_RADIANS) == 0)
1872 hpos = swe_difdeg2n(hpos, hpos2);
1873 else
1874 hpos = swe_difrad2n(hpos, hpos2);
1875 } else { /* DIFF_MIDP */
1876 if ((iflag & SEFLG_RADIANS) == 0)
1877 hpos = swe_deg_midp(hpos, hpos2);
1878 else
1879 hpos = swe_rad_midp(hpos, hpos2);
1880 }
1881 }
1882 }
1883 strcpy(spnam, se_pname);
1884 print_line(0, is_first, 0);
1885 is_first = FALSE;
1886 if (! list_hor) line_count++;
1887 if (do_orbital_elements) {
1888 orbital_elements(te, ipl, iflag, serr);
1889 continue;
1890 }
1891 if (line_count >= line_limit) {
1892 printf("****** line count %d was exceeded\n", line_limit);
1893 break;
1894 }
1895 } /* for psp */
1896 if (list_hor) {
1897 printf("\n");
1898 line_count++;
1899 }
1900 if (do_houses) {
1901 double cusp[37];
1902 double cusp_speed[37];
1903 double ascmc[10];
1904 double ascmc_speed[10];
1905 int iofs;
1906 if (toupper(ihsy) == 'G') // Gauquelin has 36 cusps
1907 nhouses = 36;
1908 iofs = nhouses + 1;
1909 iflgret = swe_houses_ex2(t,iflag, top_lat, top_long, ihsy, cusp, ascmc, cusp_speed, ascmc_speed, serr);
1910 // when swe_houses_ex() fails (e.g. with Placidus, Gauquelin, Makranski),
1911 // it always returns Porphyry cusps instead
1912 if (iflgret < 0) {
1913 char *shsy = swe_house_name(ihsy);
1914 sprintf(serr, "House method %s failed, Porphyry calculated instead", shsy);
1915 if (strcmp(serr, serr_save) != 0 ) {
1916 fputs("error: ", stdout);
1917 fputs(serr, stdout);
1918 fputs("\n", stdout);
1919 }
1920 strcpy(serr_save, serr);
1921 ihsy = 'O';
1922 nhouses = 12; // instead of 36 with 'G'
1923 iofs = nhouses + 1;
1924 }
1925 is_first = TRUE;
1926 for (ipl = 1; ipl < iofs+8; ipl++) {
1927 x[0] = cusp[ipl];
1928 if (ipl >= iofs) {
1929 x[0] = ascmc[ipl - iofs];
1930 x[3] = ascmc_speed[ipl - iofs];
1931 } else {
1932 x[3] = cusp_speed[ipl];
1933 }
1934 x[1] = 0; /* latitude */
1935 x[2] = 1.0; /* pseudo radius vector */
1936 if (ipl == iofs+2) { /* armc is already equatorial! */
1937 xequ[0] = x[0];
1938 xequ[1] = x[1];
1939 xequ[2] = x[2];
1940 } else if (strpbrk(fmt, "aADdQ") != NULL) {
1941 swe_cotrans(x, xequ, -xobl[0]);
1942 }
1943 if (strpbrk(fmt, "IiHhKk") != NULL) {
1944 double gpos[3];
1945 gpos[0] = top_long;
1946 gpos[1] = top_lat;
1947 gpos[2] = 0;
1948 swe_azalt(t, SE_ECL2HOR, gpos, datm[0], datm[1], x, xaz);
1949 }
1950 if (strpbrk(fmt, "gGj") != NULL) {
1951 hposj = swe_house_pos(armc, geopos[1], xobl[0], ihsy, x, serr);
1952 if (toupper(ihsy) == 'G')
1953 hpos = (hposj - 1) * 10;
1954 else
1955 hpos = (hposj - 1) * 30;
1956 }
1957 print_line(MODE_HOUSE, is_first, 0);
1958 is_first = FALSE;
1959 if (! list_hor) line_count++;
1960 }
1961 if (list_hor) {
1962 printf("\n");
1963 line_count++;
1964 }
1965 }
1966 if (line_count >= line_limit) {
1967 printf("****** line count %d was exceeded\n", line_limit);
1968 break;
1969 }
1970 } /* for tjd */
1971 if (*serr_warn != '\0') {
1972 printf("\nwarning: ");
1973 fputs(serr_warn,stdout);
1974 printf("\n");
1975 }
1976 } /* while 1 */
1977 /* close open files and free allocated space */
1978 end_main:
1979 if (do_set_astro_models) {
1980 printf(smod);
1981 }
1982 swe_close();
1983 return OK;
1984 }
1985
call_swe_fixstar(char * star,double te,int32 iflag,double * x,char * serr)1986 static int32 call_swe_fixstar(char *star, double te, int32 iflag, double *x, char *serr)
1987 {
1988 if (use_swe_fixstar2)
1989 return swe_fixstar2(star, te, iflag, x, serr);
1990 else
1991 return swe_fixstar(star, te, iflag, x, serr);
1992 }
1993
1994 /* This function calculates the geocentric relative distance of a planet,
1995 * where the closest position has value 1000, and remotest position has
1996 * value 0.
1997 * The value is returned as an integer. The algorithm does not allow
1998 * much higher accuracy.
1999 *
2000 * With the Moon we measure the distance relative to the maximum and minimum
2001 * found between 12000 BCE and 16000 CE.
2002 * If the distance value were given relative to the momentary osculating
2003 * ellipse, then the apogee would always have the value 1000 and the perigee
2004 * the value 0. It is certainly more interesting to know how much it is
2005 * relative to a greater time range.
2006 */
get_geocentric_relative_distance(double tjd_et,int32 ipl,int32 iflag,char * serr)2007 static int32 get_geocentric_relative_distance(double tjd_et, int32 ipl, int32 iflag, char *serr)
2008 {
2009 int32 iflagi = (iflag & (SEFLG_EPHMASK | SEFLG_HELCTR | SEFLG_BARYCTR));
2010 int32 retval;
2011 double ar = 0;
2012 double xx[6];
2013 double dmax, dmin, dtrue;
2014 if ((0) && ipl == SE_MOON) {
2015 dmax = 0.002718774; // jd = 283030.8
2016 dmin = 0.002381834; // jd = -1006731.3
2017 if ((retval = swe_calc(tjd_et, SE_MOON, iflagi | SEFLG_J2000 | SEFLG_TRUEPOS, xx, serr)) == ERR)
2018 return 0;
2019 dtrue = xx[2];
2020 } else {
2021 if (swe_orbit_max_min_true_distance(tjd_et, ipl, iflagi, &dmax, &dmin, &dtrue, serr) == ERR)
2022 return 0;
2023 }
2024 if (dmax - dmin == 0) {
2025 ar = 0;
2026 } else {
2027 ar = (1 - (dtrue - dmin) / (dmax - dmin)) * 1000.0;
2028 ar += 0.5; // rounding
2029 }
2030 return (int32) ar;
2031 }
2032
2033 /*
2034 * The string fmt contains a sequence of format specifiers;
2035 * each character in fmt creates a column, the columns are
2036 * sparated by the gap string.
2037 * Time columns tTJyY are only printed, if is_first is TRUE,
2038 * so that they are not repeated in list_hor (horizontal list) mode.
2039 * In list_hor mode, no newline is printed.
2040 */
print_line(int mode,AS_BOOL is_first,int sid_mode)2041 static int print_line(int mode, AS_BOOL is_first, int sid_mode)
2042 {
2043 char *sp, *sp2;
2044 double t2, ju2 = 0;
2045 double y_frac;
2046 double ar, sinp;
2047 double dret[20];
2048 char slon[20];
2049 char pnam[30];
2050 AS_BOOL is_house = ((mode & MODE_HOUSE) != 0);
2051 AS_BOOL is_label = ((mode & MODE_LABEL) != 0);
2052 AS_BOOL is_ayana = ((mode & MODE_AYANAMSA) != 0);
2053 int32 iflgret, dar;
2054 // build planet name column, just in case
2055 if (is_house) {
2056 if (ipl <= nhouses) {
2057 sprintf(pnam, "house %2d ", ipl);
2058 } else {
2059 sprintf(pnam, "%-15s", hs_nam[ipl - nhouses]);
2060 }
2061 } else if (diff_mode == DIFF_DIFF) {
2062 sprintf(pnam, "%.3s-%.3s", spnam, spnam2);
2063 } else if (diff_mode == DIFF_GEOHEL) {
2064 sprintf(pnam, "%.3s-%.3sHel", spnam, spnam2);
2065 } else if (diff_mode == DIFF_MIDP) {
2066 sprintf(pnam, "%.3s/%.3s", spnam, spnam2);
2067 } else {
2068 sprintf(pnam, "%-15s", spnam);
2069 }
2070 if (list_hor && strchr(fmt, 'P') == NULL) {
2071 sprintf(slon, "%.8s %s", pnam, "long.");
2072 } else {
2073 sprintf(slon, "%-14s", "long.");
2074 }
2075 for (sp = fmt; *sp != '\0'; sp++) {
2076 // if (is_house && ipl <= nhouses && strchr("bBsSrRxXuUQnNfFj+-*/=", *sp) != NULL) continue;
2077 if (is_house && strchr("bBrRxXuUQnNfFj+-*/=", *sp) != NULL) continue;
2078 if (is_ayana && strchr("bBsSrRxXuUQnNfFj+-*/=", *sp) != NULL) continue;
2079 if (sp != fmt)
2080 fputs(gap,stdout);
2081 if (sp == fmt && list_hor && !is_first && strchr("yYJtT", *sp) == NULL)
2082 fputs(gap,stdout);
2083 switch(*sp) {
2084 case 'y':
2085 if (list_hor && ! is_first) {
2086 break;
2087 }
2088 if (is_label) { printf("year"); break; }
2089 printf("%d", jyear);
2090 break;
2091 case 'Y':
2092 if (list_hor && ! is_first) {
2093 break;
2094 }
2095 if (is_label) { printf("year"); break; }
2096 t2 = swe_julday(jyear,1,1,ju2,gregflag);
2097 y_frac = (t - t2) / 365.0;
2098 printf("%.2f", jyear + y_frac);
2099 break;
2100 case 'p':
2101 if (is_label) { printf("obj.nr"); break; }
2102 if (! is_house && diff_mode == DIFF_DIFF) {
2103 printf("%d-%d", ipl, ipldiff);
2104 } else if (! is_house && diff_mode == DIFF_GEOHEL) {
2105 printf("%d-%dhel", ipl, ipldiff);
2106 } else if (! is_house && diff_mode == DIFF_MIDP) {
2107 printf("%d/%d", ipl, ipldiff);
2108 } else {
2109 printf("%d", ipl);
2110 }
2111 break;
2112 case 'P':
2113 if (is_label) { printf("%-15s", "name"); break; }
2114 if (is_house) {
2115 if (ipl <= nhouses) {
2116 printf("house %2d ", ipl);
2117 } else {
2118 printf("%-15s", hs_nam[ipl - nhouses]);
2119 }
2120 } else if (is_ayana) {
2121 // printf("Ayanamsha ");
2122 printf("Ayanamsha %s ", swe_get_ayanamsa_name(sid_mode));
2123 } else if (diff_mode == DIFF_DIFF || diff_mode == DIFF_GEOHEL) {
2124 printf("%.3s-%.3s", spnam, spnam2);
2125 } else if (diff_mode == DIFF_MIDP) {
2126 printf("%.3s/%.3s", spnam, spnam2);
2127 } else {
2128 printf("%-15s", spnam);
2129 }
2130 break;
2131 case 'J':
2132 if (list_hor && ! is_first) {
2133 break;
2134 }
2135 if (is_label) { printf("julday"); break; }
2136 y_frac = (t - floor(t)) * 100;
2137 if (floor(y_frac) != y_frac) {
2138 printf("%.5f", t);
2139 } else {
2140 printf("%.2f", t);
2141 }
2142 break;
2143 case 'T':
2144 if (list_hor && ! is_first) {
2145 break;
2146 }
2147 if (is_label) { printf("date "); break; }
2148 printf("%02d.%02d.%d", jday, jmon, jyear);
2149 if (jut != 0 || step_in_minutes || step_in_seconds ) {
2150 int h, m, s;
2151 s = (int) (jut * 3600 + 0.5);
2152 h = (int) (s / 3600.0);
2153 m = (int) ((s % 3600) / 60.0);
2154 s %= 60;
2155 printf(" %d:%02d:%02d", h, m, s);
2156 if (universal_time)
2157 printf(" UT");
2158 else
2159 printf(" TT");
2160 }
2161 break;
2162 case 't':
2163 if (list_hor && ! is_first) {
2164 break;
2165 }
2166 if (is_label) { printf("date"); break; }
2167 printf("%02d%02d%02d", jyear % 100, jmon, jday);
2168 break;
2169 case 'L':
2170 if (is_label) { printf(slon); break; }
2171 if (psp != NULL && (*psp == 'q' || *psp == 'y')) { /* delta t or time equation */
2172 printf("%# 11.7f", x[0]);
2173 printf("s");
2174 break;
2175 }
2176 fputs(dms(x[0], round_flag),stdout);
2177 break;
2178 case 'l':
2179 if (is_label) { printf(slon); break; }
2180 if (output_extra_prec)
2181 printf("%# 11.11f", x[0]);
2182 else
2183 printf("%# 11.7f", x[0]);
2184 break;
2185 case 'G':
2186 if (is_label) { printf("housPos"); break; }
2187 fputs(dms(hpos, round_flag),stdout);
2188 break;
2189 case 'g':
2190 if (is_label) { printf("housPos"); break; }
2191 printf("%# 11.7f", hpos);
2192 break;
2193 case 'j':
2194 if (is_label) { printf("houseNr"); break; }
2195 printf("%# 11.7f", hposj);
2196 break;
2197 case 'Z':
2198 if (is_label) { printf(slon); break; }
2199 fputs(dms(x[0], round_flag|BIT_ZODIAC),stdout);
2200 break;
2201 case 'S':
2202 case 's':
2203 if (*(sp+1) == 'S' || *(sp+1) == 's' || strpbrk(fmt, "XUxu") != NULL) {
2204 for (sp2 = fmt; *sp2 != '\0'; sp2++) {
2205 if (sp2 != fmt)
2206 fputs(gap,stdout);
2207 switch(*sp2) {
2208 case 'L': /* speed! */
2209 case 'Z': /* speed! */
2210 if (is_label) { printf("lon/day"); break; }
2211 fputs(dms(x[3], round_flag),stdout);
2212 break;
2213 case 'l': /* speed! */
2214 if (is_label) { printf("lon/day"); break; }
2215 if (output_extra_prec)
2216 printf("%# 11.9f", x[3]);
2217 else
2218 printf("%# 11.7f", x[3]);
2219 break;
2220 case 'B': /* speed! */
2221 if (is_label) { printf("lat/day"); break; }
2222 fputs(dms(x[4], round_flag),stdout);
2223 break;
2224 case 'b': /* speed! */
2225 if (is_label) { printf("lat/day"); break; }
2226 if (output_extra_prec)
2227 printf("%# 11.9f", x[4]);
2228 else
2229 printf("%# 11.7f", x[4]);
2230 break;
2231 case 'A': /* speed! */
2232 if (is_label) { printf("RA/day"); break; }
2233 fputs(dms(xequ[3]/15, round_flag|SEFLG_EQUATORIAL),stdout);
2234 break;
2235 case 'a': /* speed! */
2236 if (is_label) { printf("RA/day"); break; }
2237 if (output_extra_prec)
2238 printf("%# 11.9f", xequ[3]);
2239 else
2240 printf("%# 11.7f", xequ[3]);
2241 break;
2242 case 'D': /* speed! */
2243 if (is_label) { printf("dcl/day"); break; }
2244 fputs(dms(xequ[4], round_flag),stdout);
2245 break;
2246 case 'd': /* speed! */
2247 if (is_label) { printf("dcl/day"); break; }
2248 if (output_extra_prec)
2249 printf("%# 11.9f", xequ[4]);
2250 else
2251 printf("%# 11.7f", xequ[4]);
2252 break;
2253 case 'R': /* speed! */
2254 case 'r': /* speed! */
2255 if (is_label) { printf("AU/day"); break; }
2256 if (output_extra_prec)
2257 printf("%# 16.14f", x[5]);
2258 else
2259 printf("%# 14.9f", x[5]);
2260 break;
2261 case 'U': /* speed! */
2262 case 'X': /* speed! */
2263 if (is_label) {
2264 fputs("speed_0", stdout);
2265 fputs(gap, stdout);
2266 fputs("speed_1", stdout);
2267 fputs(gap, stdout);
2268 fputs("speed_2", stdout);
2269 break;
2270 }
2271 if (*sp =='U')
2272 ar = sqrt(square_sum(xcart));
2273 else
2274 ar = 1;
2275 printf("%# 14.9f", xcart[3]/ar);
2276 fputs(gap,stdout);
2277 printf("%# 14.9f", xcart[4]/ar);
2278 fputs(gap,stdout);
2279 printf("%# 14.9f", xcart[5]/ar);
2280 break;
2281 case 'u': /* speed! */
2282 case 'x': /* speed! */
2283 if (is_label) {
2284 fputs("speed_0", stdout);
2285 fputs(gap, stdout);
2286 fputs("speed_1", stdout);
2287 fputs(gap, stdout);
2288 fputs("speed_2", stdout);
2289 break;
2290 }
2291 if (*sp =='u')
2292 ar = sqrt(square_sum(xcartq));
2293 else
2294 ar = 1;
2295 printf("%# 14.9f", xcartq[3]/ar);
2296 fputs(gap,stdout);
2297 printf("%# 14.9f", xcartq[4]/ar);
2298 fputs(gap,stdout);
2299 printf("%# 14.9f", xcartq[5]/ar);
2300 break;
2301 default:
2302 break;
2303 }
2304 }
2305 if (*(sp+1) == 'S' || *(sp+1) == 's')
2306 sp++;
2307 } else if (*sp == 'S') {
2308 int flag = round_flag;
2309 if (is_house) flag |= BIT_ALLOW_361; // speed of houses can be > 360
2310 if (is_label) { printf("deg/day"); break; }
2311 fputs(dms(x[3], flag),stdout);
2312 } else {
2313 if (is_label) { printf("deg/day"); break; }
2314 if (output_extra_prec)
2315 printf("%# 11.17f", x[3]);
2316 else
2317 printf("%# 11.7f", x[3]);
2318 }
2319 break;
2320 case 'B':
2321 if (is_label) { printf("lat. "); break; }
2322 if (*psp == 'q') { /* delta t */
2323 printf("%# 11.7f", x[1]);
2324 printf("h");
2325 break;
2326 }
2327 fputs(dms(x[1], round_flag),stdout);
2328 break;
2329 case 'b':
2330 if (is_label) { printf("lat. "); break; }
2331 if (output_extra_prec)
2332 printf("%# 11.11f", x[1]);
2333 else
2334 printf("%# 11.7f", x[1]);
2335 break;
2336 case 'A': /* right ascension */
2337 if (is_label) { printf("RA "); break; }
2338 fputs(dms(xequ[0]/15, round_flag|SEFLG_EQUATORIAL),stdout);
2339 break;
2340 case 'a': /* right ascension */
2341 if (is_label) { printf("RA "); break; }
2342 if (output_extra_prec)
2343 printf("%# 11.11f", xequ[0]);
2344 else
2345 printf("%# 11.7f", xequ[0]);
2346 break;
2347 case 'D': /* declination */
2348 if (is_label) { printf("decl "); break; }
2349 fputs(dms(xequ[1], round_flag),stdout);
2350 break;
2351 case 'd': /* declination */
2352 if (is_label) { printf("decl "); break; }
2353 if (output_extra_prec)
2354 printf("%# 11.11f", xequ[1]);
2355 else
2356 printf("%# 11.7f", xequ[1]);
2357 break;
2358 case 'I': /* azimuth */
2359 if (is_label) { printf("azimuth"); break; }
2360 fputs(dms(xaz[0], round_flag),stdout);
2361 break;
2362 case 'i': /* azimuth */
2363 if (is_label) { printf("azimuth"); break; }
2364 printf("%# 11.7f", xaz[0]);
2365 break;
2366 case 'H': /* height */
2367 if (is_label) { printf("height"); break; }
2368 fputs(dms(xaz[1], round_flag),stdout);
2369 break;
2370 case 'h': /* height */
2371 if (is_label) { printf("height"); break; }
2372 printf("%# 11.7f", xaz[1]);
2373 break;
2374 case 'K': /* height (apparent) */
2375 if (is_label) { printf("hgtApp"); break; }
2376 fputs(dms(xaz[2], round_flag),stdout);
2377 break;
2378 case 'k': /* height (apparent) */
2379 if (is_label) { printf("hgtApp"); break; }
2380 printf("%# 11.7f", xaz[2]);
2381 break;
2382 case 'R':
2383 if (is_label) { printf("distAU "); break; }
2384 if (output_extra_prec)
2385 printf("%# 16.14f", x[2]);
2386 else
2387 printf("%# 14.9f", x[2]);
2388 break;
2389 case 'W':
2390 if (is_label) { printf("distLY "); break; }
2391 printf("%# 14.9f", x[2] * SE_AUNIT_TO_LIGHTYEAR);
2392 break;
2393 case 'w':
2394 if (is_label) { printf("distkm "); break; }
2395 printf("%# 14.9f", x[2] * SE_AUNIT_TO_KM);
2396 break;
2397 case 'r':
2398 if (is_label) { printf("dist"); break; }
2399 if ( ipl == SE_MOON ) { /* for moon print parallax */
2400 /* geocentric horizontal parallax: */
2401 if ((0)) {
2402 sinp = 8.794 / x[2]; /* in seconds of arc */
2403 ar = sinp * (1 + sinp * sinp * 3.917402e-12);
2404 /* the factor is 1 / (3600^2 * (180/pi)^2 * 6) */
2405 printf("%# 13.5f\" %# 13.5f'", ar, ar/60.0);
2406 }
2407 swe_pheno(te, ipl, iflag, dret, serr);
2408 printf("%# 13.5f\"", dret[5] * 3600);
2409 } else {
2410 printf("%# 14.9f", x[2]);
2411 }
2412 break;
2413 case 'q':
2414 if (is_label) { printf("reldist"); break; }
2415 dar = get_geocentric_relative_distance(te, ipl, iflag, serr);
2416 printf("% 5d", dar);
2417 break;
2418 case 'U':
2419 case 'X':
2420 if (*sp =='U')
2421 ar = sqrt(square_sum(xcart));
2422 else
2423 ar = 1;
2424 printf("%# 14.9f", xcart[0]/ar);
2425 fputs(gap,stdout);
2426 printf("%# 14.9f", xcart[1]/ar);
2427 fputs(gap,stdout);
2428 printf("%# 14.9f", xcart[2]/ar);
2429 break;
2430 case 'u':
2431 case 'x':
2432 if (is_label) {
2433 fputs("x0", stdout);
2434 fputs(gap, stdout);
2435 fputs("x1", stdout);
2436 fputs(gap, stdout);
2437 fputs("x2", stdout);
2438 break;
2439 }
2440 if (*sp =='u')
2441 ar = sqrt(square_sum(xcartq));
2442 else
2443 ar = 1;
2444 if (output_extra_prec) {
2445 printf("%# .17f", xcartq[0]/ar);
2446 fputs(gap,stdout);
2447 printf("%# .17f", xcartq[1]/ar);
2448 fputs(gap,stdout);
2449 printf("%# .17f", xcartq[2]/ar);
2450 } else {
2451 printf("%# 14.9f", xcartq[0]/ar);
2452 fputs(gap,stdout);
2453 printf("%# 14.9f", xcartq[1]/ar);
2454 fputs(gap,stdout);
2455 printf("%# 14.9f", xcartq[2]/ar);
2456 }
2457 break;
2458 case 'Q':
2459 if (is_label) { printf("Q"); break; }
2460 printf("%-15s", spnam);
2461 fputs(dms(x[0], round_flag),stdout);
2462 fputs(dms(x[1], round_flag),stdout);
2463 printf(" %# 14.9f", x[2]);
2464 fputs(dms(x[3], round_flag),stdout);
2465 fputs(dms(x[4], round_flag),stdout);
2466 printf(" %# 14.9f\n", x[5]);
2467 printf(" %s", dms(xequ[0], round_flag));
2468 fputs(dms(xequ[1], round_flag),stdout);
2469 printf(" %s", dms(xequ[3], round_flag));
2470 fputs(dms(xequ[4], round_flag),stdout);
2471 break;
2472 case 'N':
2473 case 'n': {
2474 double xasc[6], xdsc[6];
2475 int imeth = (*sp == tolower(*sp))?SE_NODBIT_MEAN:SE_NODBIT_OSCU;
2476 iflgret = swe_nod_aps(te, ipl, iflag, imeth, xasc, xdsc, NULL, NULL, serr);
2477 if (iflgret >= 0 && (ipl <= SE_NEPTUNE || *sp == 'N') ) {
2478 if (is_label) {
2479 fputs("nodAsc", stdout);
2480 fputs(gap, stdout);
2481 fputs("nodDesc", stdout);
2482 break;
2483 }
2484 if (use_dms)
2485 fputs(dms(xasc[0], round_flag|BIT_ZODIAC),stdout);
2486 else
2487 printf("%# 11.7f", xasc[0]);
2488 fputs(gap,stdout);
2489 if (use_dms)
2490 fputs(dms(xdsc[0], round_flag|BIT_ZODIAC),stdout);
2491 else
2492 printf("%# 11.7f", xdsc[0]);
2493 }
2494 };
2495 break;
2496 case 'F':
2497 case 'f':
2498 if (! is_house) {
2499 double xfoc[6], xaph[6], xper[6];
2500 int imeth = (*sp == tolower(*sp))?SE_NODBIT_MEAN:SE_NODBIT_OSCU;
2501 // fprintf(stderr, "c=%c\n", *sp);
2502 iflgret = swe_nod_aps(te, ipl, iflag, imeth, NULL, NULL, xper, xaph, serr);
2503 if (iflgret >= 0 && (ipl <= SE_NEPTUNE || *sp == 'F') ) {
2504 if (is_label) {
2505 fputs("peri", stdout);
2506 fputs(gap, stdout);
2507 fputs("apo", stdout);
2508 fputs(gap, stdout);
2509 fputs("focus", stdout);
2510 break;
2511 }
2512 printf("%# 11.7f", xper[0]);
2513 fputs(gap,stdout);
2514 printf("%# 11.7f", xaph[0]);
2515 }
2516 imeth |= SE_NODBIT_FOPOINT;
2517 iflgret = swe_nod_aps(te, ipl, iflag, imeth, NULL, NULL, xper, xfoc, serr);
2518 if (iflgret >= 0 && (ipl <= SE_NEPTUNE || *sp == 'F') ) {
2519 fputs(gap,stdout);
2520 printf("%# 11.7f", xfoc[0]);
2521 }
2522 };
2523 break;
2524 case '+':
2525 if (is_house) break;
2526 if (is_label) { printf("phase"); break; }
2527 fputs(dms(attr[0], round_flag),stdout);
2528 break;
2529 case '-':
2530 if (is_label) { printf("phase"); break; }
2531 if (is_house) break;
2532 printf(" %# 14.9f", attr[1]);
2533 break;
2534 case '*':
2535 if (is_label) { printf("elong"); break; }
2536 if (is_house) break;
2537 fputs(dms(attr[2], round_flag),stdout);
2538 break;
2539 case '/':
2540 if (is_label) { printf("diamet"); break; }
2541 if (is_house) break;
2542 fputs(dms(attr[3], round_flag),stdout);
2543 break;
2544 case '=':
2545 if (is_label) { printf("magn"); break; }
2546 if (is_house) break;
2547 printf(" %# 6.3fm", attr[4]);
2548 break;
2549 case 'V': /* human design gates */
2550 case 'v': {
2551 double xhds;
2552 int igate, iline, ihex;
2553 static int hexa[64] = {1, 43, 14, 34, 9, 5, 26, 11, 10, 58, 38, 54, 61, 60, 41, 19, 13, 49, 30, 55, 37, 63, 22, 36, 25, 17, 21, 51, 42, 3, 27, 24, 2, 23, 8, 20, 16, 35, 45, 12, 15, 52, 39, 53, 62, 56, 31, 33, 7, 4, 29, 59, 40, 64, 47, 6, 46, 18, 48, 57, 32, 50, 28, 44};
2554 if (is_label) { printf("hds"); break; }
2555 if (is_house) break;
2556 xhds = swe_degnorm(x[0] - 223.25);
2557 ihex = (int) floor(xhds / 5.625);
2558 iline = ((int) (floor(xhds / 0.9375))) % 6 + 1 ;
2559 igate = hexa[ihex];
2560 printf("%2d.%d", igate, iline);
2561 if (*sp == 'V')
2562 printf(" %2d%%", swe_d2l(100 * fmod(xhds / 0.9375, 1)));
2563 break;
2564 }
2565 case 'm': { // Meridian distance
2566 if (is_label) { printf("MD "); break; }
2567 double md = swe_difdeg2n(xequ[0], armc);
2568 if (md < 0) md = -md;
2569 if (output_extra_prec)
2570 printf("%# 11.11f", md);
2571 else
2572 printf("%# 11.7f", md);
2573 break;
2574 }
2575 case 'z': { // Zenith distance
2576 if (is_label) { printf("ZD "); break; }
2577 swe_azalt(tut, SE_EQU2HOR, geopos, datm[0], datm[1], xequ, xaz);
2578 double zd = 90 - xaz[1];
2579 if (output_extra_prec)
2580 printf("%# 11.11f", zd);
2581 else
2582 printf("%# 11.7f", zd);
2583 break;
2584 }
2585 } /* switch */
2586 } /* for sp */
2587 if (! list_hor)
2588 printf("\n");
2589 return OK;
2590 }
2591
dms(double xv,int32 iflg)2592 static char *dms(double xv, int32 iflg)
2593 {
2594 int izod;
2595 int32 k, kdeg, kmin, ksec;
2596 char *c = ODEGREE_STRING;
2597 char *sp, s1[50];
2598 static char s[50];
2599 int sgn;
2600 #if MSDOS
2601 if (_isnan(xv))
2602 return "nan";
2603 #else
2604 if (isnan(xv))
2605 return "nan";
2606 #endif
2607 if (xv >= 360 && !(iflg & BIT_ALLOW_361))
2608 xv = 0;
2609 *s = '\0';
2610 if (iflg & SEFLG_EQUATORIAL)
2611 c = "h";
2612 if (xv < 0) {
2613 xv = -xv;
2614 sgn = -1;
2615 } else {
2616 sgn = 1;
2617 }
2618 if (iflg & BIT_ROUND_MIN) {
2619 if (!(iflg & BIT_ALLOW_361))
2620 xv = swe_degnorm(xv + 0.5/60);
2621 } else if (iflg & BIT_ROUND_SEC) {
2622 if (!(iflg & BIT_ALLOW_361))
2623 xv = swe_degnorm(xv + 0.5/3600);
2624 } else {
2625 /* rounding 0.9999999999 to 1 */
2626 if (output_extra_prec)
2627 xv += (xv < 0 ? -1 : 1 ) * 0.000000005 / 3600.0;
2628 else
2629 xv += (xv < 0 ? -1 : 1 ) * 0.00005 / 3600.0;
2630 }
2631 if (iflg & BIT_ZODIAC) {
2632 izod = (int) (xv / 30);
2633 if (izod == 12) izod = 0;
2634 xv = fmod(xv, 30);
2635 kdeg = (int32) xv;
2636 sprintf(s, "%2d %s ", kdeg, zod_nam[izod]);
2637 } else {
2638 kdeg = (int32) xv;
2639 sprintf(s, " %3d%s", kdeg, c);
2640 }
2641 xv -= kdeg;
2642 xv *= 60;
2643 kmin = (int32) xv;
2644 if ((iflg & BIT_ZODIAC) && (iflg & BIT_ROUND_MIN)) {
2645 sprintf(s1, "%2d", kmin);
2646 } else {
2647 sprintf(s1, "%2d'", kmin);
2648 }
2649 strcat(s, s1);
2650 if (iflg & BIT_ROUND_MIN)
2651 goto return_dms;
2652 xv -= kmin;
2653 xv *= 60;
2654 ksec = (int32) xv;
2655 if (iflg & BIT_ROUND_SEC) {
2656 sprintf(s1, "%2d\"", ksec);
2657 } else {
2658 sprintf(s1, "%2d", ksec);
2659 }
2660 strcat(s, s1);
2661 if (iflg & BIT_ROUND_SEC)
2662 goto return_dms;
2663 xv -= ksec;
2664 if (output_extra_prec) {
2665 k = (int32) (xv * 100000000);
2666 sprintf(s1, ".%08d", k);
2667 } else {
2668 k = (int32) (xv * 10000);
2669 sprintf(s1, ".%04d", k);
2670 }
2671 strcat(s, s1);
2672 return_dms:;
2673 if (sgn < 0) {
2674 sp = strpbrk(s, "0123456789");
2675 *(sp-1) = '-';
2676 }
2677 if (iflg & BIT_LZEROES) {
2678 while ((sp = strchr(s+2, ' ')) != NULL) *sp = '0';
2679 }
2680 return(s);
2681 }
2682
letter_to_ipl(int letter)2683 static int letter_to_ipl(int letter)
2684 {
2685 if (letter >= '0' && letter <= '9')
2686 return letter - '0' + SE_SUN;
2687 if (letter >= 'A' && letter <= 'I')
2688 return letter - 'A' + SE_MEAN_APOG;
2689 if (letter >= 'J' && letter <= 'Z')
2690 return letter - 'J' + SE_CUPIDO;
2691 switch (letter) {
2692 case 'm': return SE_MEAN_NODE;
2693 case 'c': return SE_INTP_APOG;
2694 case 'g': return SE_INTP_PERG;
2695 case 'n':
2696 case 'o': return SE_ECL_NUT;
2697 case 't': return SE_TRUE_NODE;
2698 case 'f': return SE_FIXSTAR;
2699 case 'w': return SE_WALDEMATH;
2700 case 'e': /* swetest: a line of labels */
2701 case 'q': /* swetest: delta t */
2702 case 'y': /* swetest: time equation */
2703 case 'x': /* swetest: sidereal time */
2704 case 'b': /* swetest: ayanamsha */
2705 case 's': /* swetest: an asteroid, with number given in -xs[number] */
2706 case 'v': /* swetest: a planetary moon, with number given in -xv[number] */
2707 case 'z': /* swetest: a fictitious body, number given in -xz[number] */
2708 case 'd': /* swetest: default (main) factors 0123456789mtABC */
2709 case 'p': /* swetest: main factors ('d') plus main asteroids DEFGHI */
2710 case 'h': /* swetest: fictitious factors JKLMNOPQRSTUVWXYZw */
2711 case 'a': /* swetest: all factors, like 'p'+'h' */
2712 return -1;
2713 }
2714 return -2;
2715 }
2716
ut_to_lmt_lat(double t_ut,double * geopos,double * t_ret,char * serr)2717 static int32 ut_to_lmt_lat(double t_ut, double *geopos, double *t_ret, char *serr)
2718 {
2719 int32 iflgret = OK;
2720 if (time_flag & (BIT_TIME_LMT | BIT_TIME_LAT)) {
2721 t_ut += geopos[0] / 360.0;
2722 if (time_flag & BIT_TIME_LAT) {
2723 iflgret = swe_lmt_to_lat(t_ut, geopos[0], &t_ut, serr);
2724 }
2725 }
2726 *t_ret = t_ut;
2727 return iflgret;
2728 }
2729
orbital_elements(double tjd_et,int32 ipl,int32 iflag,char * serr)2730 static int32 orbital_elements(double tjd_et, int32 ipl, int32 iflag, char *serr)
2731 {
2732 int32 retval;
2733 double dret[20], jut;
2734 int32 jyear, jmon, jday;
2735 char sdateperi[20];
2736 retval = swe_get_orbital_elements(tjd_et, ipl, iflag, dret, serr);
2737 if (retval == ERR) {
2738 printf("%s\n", serr);
2739 return ERR;
2740 } else {
2741 swe_revjul(dret[14], gregflag, &jyear, &jmon, &jday, &jut);
2742 sprintf(sdateperi, "%2d.%02d.%04d,%s", jday, jmon, jyear, hms(jut,BIT_LZEROES));
2743 printf("semiaxis \t%f\neccentricity \t%f\ninclination \t%f\nasc. node \t%f\narg. pericenter \t%f\npericenter \t%f\n", dret[0], dret[1], dret[2], dret[3], dret[4], dret[5]);
2744 printf("mean longitude \t%f\nmean anomaly \t%f\necc. anomaly \t%f\ntrue anomaly \t%f\n", dret[9], dret[6], dret[8], dret[7]);
2745 printf("time pericenter \t%f %s\ndist. pericenter \t%f\ndist. apocenter \t%f\n", dret[14], sdateperi, dret[15], dret[16]);
2746 printf("mean daily motion\t%f\nsid. period (y) \t%f\ntrop. period (y) \t%f\nsynodic cycle (d)\t%f\n", dret[11], dret[10], dret[12], dret[13]);
2747 }
2748 return OK;
2749 }
2750
insert_gap_string_for_tabs(char * sout,char * gap)2751 static void insert_gap_string_for_tabs(char *sout, char *gap)
2752 {
2753 char *sp;
2754 char s[LEN_SOUT];
2755 if (!have_gap_parameter)
2756 return;
2757 if (*gap == '\t')
2758 return;
2759 while((sp = strchr(sout, '\t')) != NULL && strlen(sout) + strlen(gap) < LEN_SOUT) {
2760 strcpy(s, sp + 1);
2761 strcpy(sp, gap);
2762 strcat(sp, s);
2763 }
2764 }
2765
print_rise_set_line(double trise,double tset,double * geopos,char * serr)2766 static int32 print_rise_set_line(double trise, double tset, double *geopos, char *serr) {
2767 double t0;
2768 int retc = OK;
2769 *sout = '\0';
2770 if (trise != 0) retc = ut_to_lmt_lat(trise, geopos, &(trise), serr);
2771 if (tset != 0) retc = ut_to_lmt_lat(tset, geopos, &(tset), serr);
2772 strcpy(sout, "rise ");
2773 if (have_gap_parameter) sprintf(sout + strlen(sout), "\t");
2774 if (trise == 0) {
2775 strcat(sout, " -\t - ");
2776 } else {
2777 swe_revjul(trise, gregflag, &jyear, &jmon, &jday, &jut);
2778 sprintf(sout + strlen(sout), "%2d.%02d.%04d\t%s ", jday, jmon, jyear, hms(jut,BIT_LZEROES));
2779 }
2780 if (have_gap_parameter) sprintf(sout + strlen(sout), "\t");
2781 strcat(sout, "set ");
2782 if (have_gap_parameter) sprintf(sout + strlen(sout), "\t");
2783 if (tset == 0 ) {
2784 strcat(sout, " -\t - ");
2785 } else {
2786 swe_revjul(tset, gregflag, &jyear, &jmon, &jday, &jut);
2787 sprintf(sout + strlen(sout), "%2d.%02d.%04d\t%s ", jday, jmon, jyear, hms(jut,BIT_LZEROES));
2788 }
2789 if (trise != 0 && tset != 0) {
2790 if (have_gap_parameter) sprintf(sout + strlen(sout), "\t");
2791 sprintf(sout + strlen(sout), "dt =");
2792 if (have_gap_parameter) sprintf(sout + strlen(sout), "\t");
2793 t0 = (tset - trise) * 24;
2794 sprintf(sout + strlen(sout), "%s", hms(t0, BIT_LZEROES));
2795 }
2796 strcat(sout, "\n");
2797 if (have_gap_parameter) insert_gap_string_for_tabs(sout, gap);
2798 do_printf(sout);
2799 return retc;
2800 }
2801
call_rise_set(double t_ut,int32 ipl,char * star,int32 whicheph,double * geopos,char * serr)2802 static int32 call_rise_set(double t_ut, int32 ipl, char *star, int32 whicheph, double *geopos, char *serr)
2803 {
2804 int ii, rval, loop_count;
2805 int32 rsmi = 0;
2806 double dayfrac = 0.0001;
2807 double tret[10], trise, tset, tnext, tret1sv = 0;
2808 AS_BOOL do_rise, do_set;
2809 AS_BOOL last_was_empty = FALSE;
2810 int32 retc = OK;
2811 int rsmior = 0;
2812 if (norefrac) rsmior |= SE_BIT_NO_REFRACTION;
2813 if (disccenter) rsmior |= SE_BIT_DISC_CENTER;
2814 if (discbottom) rsmior |= SE_BIT_DISC_BOTTOM;
2815 if (hindu) rsmior |= SE_BIT_HINDU_RISING;
2816 if (fabs(geopos[1]) < 60 && ipl >= SE_SUN && ipl <= SE_PLUTO)
2817 dayfrac = 0.01;
2818 swe_set_topo(geopos[0], geopos[1], geopos[2]);
2819 // "geo. long 8.000000, lat 47.000000, alt 0.000000"
2820 if (with_header)
2821 printf("\ngeo. long %f, lat %f, alt %f", geopos[0], geopos[1], geopos[2]);
2822 do_printf("\n");
2823 tnext = t_ut;
2824 // the code is designed for looping with -nxxx over many days, during which
2825 // the object might become circumpolar, or never rise at all.
2826 while (special_event == SP_RISE_SET && tnext < t_ut + nstep) {
2827 // the following 'if' avoids unnecessary calculations for circumpolar
2828 // objects. even without it, the output would be correct, but
2829 // could be considerably slower.
2830 if (last_was_empty && (star == NULL || *star == '\0')) {
2831 rval = swe_calc_ut(tnext, ipl, whicheph|SEFLG_EQUATORIAL, tret, serr);
2832 if (rval >= 0 ) {
2833 double edist = geopos[1] + tret[1];
2834 double edist2 = geopos[1] - tret[1];
2835 if ((edist - 2 > 90 || edist + 2 < -90)
2836 || (edist2 - 2 > 90 || edist2 + 2 < -90)) {
2837 tnext += 1;
2838 continue;
2839 }
2840 }
2841 }
2842 /* rising */
2843 rsmi = SE_CALC_RISE | rsmior;
2844 rval= swe_rise_trans(tnext, ipl, star, whicheph, rsmi, geopos, datm[0], datm[1], &trise, serr);
2845 if (rval == ERR) {
2846 do_printf(serr);
2847 exit(0);
2848 }
2849 do_rise = (rval == OK);
2850 /* setting */
2851 rsmi = SE_CALC_SET | rsmior;
2852 do_set = FALSE;
2853 loop_count = 0;
2854 //tnext = trise; // dieter 14-feb-17
2855 while (! do_set && loop_count < 2) {
2856 rval = swe_rise_trans(tnext, ipl, star, whicheph, rsmi, geopos, datm[0], datm[1], &tset, serr);
2857 if (rval == ERR) {
2858 do_printf(serr);
2859 exit(0);
2860 }
2861 do_set = (rval == OK);
2862 if (!do_set && do_rise ) {
2863 tnext = trise;
2864 }
2865 loop_count++;
2866 }
2867 if (do_rise && do_set && trise > tset) {
2868 do_rise = FALSE; // ignore rises happening before setting
2869 trise = 0; // we hope that exact time 0 never happens, is highly unlikely.
2870 }
2871 if (do_rise && do_set) {
2872 rval = print_rise_set_line(trise, tset, geopos, serr);
2873 last_was_empty = FALSE;
2874 tnext = tset + dayfrac;
2875 } else if (do_rise && !do_set) {
2876 rval = print_rise_set_line(trise, 0, geopos, serr);
2877 last_was_empty = FALSE;
2878 tnext = trise + dayfrac;
2879 } else if (do_set && ! do_rise) {
2880 tnext = tset + dayfrac;
2881 rval = print_rise_set_line(0, tset, geopos, serr);
2882 last_was_empty = FALSE;
2883 } else { // neither rise nor set
2884 // for sequences of days without rise or set, the line '- -' is printed only once.
2885 if (! last_was_empty) rval = print_rise_set_line(0, 0, geopos, serr);
2886 tnext += 1;
2887 last_was_empty = TRUE;
2888 }
2889 if (rval == ERR) {
2890 do_printf(serr);
2891 exit(0);
2892 }
2893 if (nstep == 1) break;
2894 }
2895 /* swetest -metr
2896 * calculate and print transits over meridian (midheaven and lower
2897 * midheaven */
2898 if (special_event == SP_MERIDIAN_TRANSIT) {
2899 /* loop over days */
2900 for (ii = 0; ii < nstep; ii++, t_ut = tret1sv + 0.001) {
2901 /* transit over midheaven */
2902 if (swe_rise_trans(t_ut, ipl, star, whicheph, SE_CALC_MTRANSIT, geopos, datm[0], datm[1], &(tret[0]), serr) != OK) {
2903 do_printf(serr);
2904 return ERR;
2905 }
2906 /* transit over lower midheaven */
2907 if (swe_rise_trans(t_ut, ipl, star, whicheph, SE_CALC_ITRANSIT, geopos, datm[0], datm[1], &(tret[1]), serr) != OK) {
2908 do_printf(serr);
2909 return ERR;
2910 }
2911 tret1sv = tret[1];
2912 if (time_flag & (BIT_TIME_LMT | BIT_TIME_LAT)) {
2913 retc = ut_to_lmt_lat(tret[0], geopos, &(tret[0]), serr);
2914 retc = ut_to_lmt_lat(tret[1], geopos, &(tret[1]), serr);
2915 }
2916 strcpy(sout, "mtransit ");
2917 if (have_gap_parameter) sprintf(sout + strlen(sout), "\t");
2918 if (tret[0] == 0 || tret[0] > tret[1]) strcat(sout, " -\t - ");
2919 else {
2920 swe_revjul(tret[0], gregflag, &jyear, &jmon, &jday, &jut);
2921 sprintf(sout + strlen(sout), "%2d.%02d.%04d\t%s ", jday, jmon, jyear, hms(jut,BIT_LZEROES));
2922 }
2923 if (have_gap_parameter) sprintf(sout + strlen(sout), "\t");
2924 strcat(sout, "itransit ");
2925 if (have_gap_parameter) sprintf(sout + strlen(sout), "\t");
2926 if (tret[1] == 0) strcat(sout, " -\t - \n");
2927 else {
2928 swe_revjul(tret[1], gregflag, &jyear, &jmon, &jday, &jut);
2929 sprintf(sout + strlen(sout), "%2d.%02d.%04d\t%s\n", jday, jmon, jyear, hms(jut,BIT_LZEROES));
2930 }
2931 if (have_gap_parameter) insert_gap_string_for_tabs(sout, gap);
2932 do_printf(sout);
2933 }
2934 }
2935 return retc;
2936 }
2937
get_gregjul(int gregflag,int year)2938 static char* get_gregjul(int gregflag, int year)
2939 {
2940 if (gregflag == SE_JUL_CAL) return "jul";
2941 if (year < 1700) return "greg";
2942 return "";
2943 }
2944
2945 // print lon and lat string in minute precision
format_lon_lat(char * slon,char * slat,double lon,double lat)2946 static void format_lon_lat(char *slon, char *slat, double lon, double lat)
2947 {
2948 int roundflag, ideg, imin, isec, isgn;
2949 double dsecfr;
2950 char c;
2951 roundflag = SE_SPLIT_DEG_ROUND_SEC;
2952 swe_split_deg(lon, roundflag, &ideg, &imin, &isec, &dsecfr, &isgn);
2953 c = (lon < 0) ? 'w' : 'e';
2954 sprintf(slon, "%d%c%02d%02d", abs(ideg), c, imin, isec);
2955 swe_split_deg(lat, roundflag, &ideg, &imin, &isec, &dsecfr, &isgn);
2956 c = (lat < 0) ? 's' : 'n';
2957 sprintf(slat, "%d%c%02d%02d", abs(ideg), c, imin, isec);
2958 }
2959
call_lunar_eclipse(double t_ut,int32 whicheph,int32 special_mode,double * geopos,char * serr)2960 static int32 call_lunar_eclipse(double t_ut, int32 whicheph, int32 special_mode, double *geopos, char *serr)
2961 {
2962 int i, ii, retc = OK, eclflag, ecl_type = 0;
2963 int rval, ihou, imin, isec, isgn;
2964 double dfrc, attr[30], dt, xx[6], geopos_max[3];
2965 char s1[AS_MAXCH], s2[AS_MAXCH], sout_short[AS_MAXCH], sfmt[AS_MAXCH], *styp = "none", *sgj;
2966 char slon[8], slat[8], saros[20];
2967 geopos_max[0] = 0; geopos_max[1] = 0;
2968 /* no selective eclipse type set, set all */
2969 if (with_chart_link) do_printf("<pre>");
2970 if ((search_flag & SE_ECL_ALLTYPES_LUNAR) == 0)
2971 search_flag |= SE_ECL_ALLTYPES_LUNAR;
2972 // "geo. long 8.000000, lat 47.000000, alt 0.000000"
2973 if (special_mode & SP_MODE_LOCAL) {
2974 if (with_header)
2975 printf("\ngeo. long %f, lat %f, alt %f", geopos[0], geopos[1], geopos[2]);
2976 }
2977 do_printf("\n");
2978 for (ii = 0; ii < nstep; ii++, t_ut += direction) {
2979 *sout = '\0';
2980 /* swetest -lunecl -how
2981 * type of lunar eclipse and percentage for a given time: */
2982 if (special_mode & SP_MODE_HOW) {
2983 if ((eclflag = swe_lun_eclipse_how(t_ut, whicheph, geopos, attr, serr)) ==
2984 ERR) {
2985 do_printf(serr);
2986 return ERR;
2987 } else {
2988 if (eclflag & SE_ECL_TOTAL) {
2989 ecl_type = ECL_LUN_TOTAL;
2990 strcpy(sfmt, "total lunar eclipse: %f o/o \n");
2991 } else if (eclflag & SE_ECL_PARTIAL) {
2992 ecl_type = ECL_LUN_PARTIAL;
2993 strcpy(sfmt, "partial lunar eclipse: %f o/o \n");
2994 } else if (eclflag & SE_ECL_PENUMBRAL) {
2995 ecl_type = ECL_LUN_PENUMBRAL;
2996 strcpy(sfmt, "penumbral lunar eclipse: %f o/o \n");
2997 } else {
2998 strcpy(sfmt, "no lunar eclipse \n");
2999 }
3000 strcpy(sout, sfmt);
3001 if (strchr(sfmt, '%') != NULL) {
3002 sprintf(sout, sfmt, attr[0]);
3003 }
3004 do_printf(sout);
3005 }
3006 continue;
3007 }
3008 /* swetest -lunecl
3009 * find next lunar eclipse: */
3010 /* locally visible lunar eclipse */
3011 if (special_mode & SP_MODE_LOCAL) {
3012 if ((eclflag = swe_lun_eclipse_when_loc(t_ut, whicheph, geopos,
3013 tret, attr, direction_flag, serr)) == ERR) {
3014 do_printf(serr);
3015 return ERR;
3016 }
3017 if (time_flag & (BIT_TIME_LMT | BIT_TIME_LAT)) {
3018 for (i = 0; i < 10; i++) {
3019 if (tret[i] != 0) {
3020 retc = ut_to_lmt_lat(tret[i], geopos, &(tret[i]), serr);
3021 if (retc == ERR) {
3022 do_printf(serr);
3023 return ERR;
3024 }
3025 }
3026 }
3027 }
3028 t_ut = tret[0];
3029 if ((eclflag & SE_ECL_TOTAL)) {
3030 strcpy(sout, "total ");
3031 ecl_type = ECL_LUN_TOTAL;
3032 }
3033 if ((eclflag & SE_ECL_PENUMBRAL)) {
3034 strcpy(sout, "penumb. ");
3035 ecl_type = ECL_LUN_PENUMBRAL;
3036 }
3037 if ((eclflag & SE_ECL_PARTIAL)) {
3038 strcpy(sout, "partial ");
3039 ecl_type = ECL_LUN_PARTIAL;
3040 }
3041 strcat(sout, "lunar eclipse\t");
3042 swe_revjul(t_ut, gregflag, &jyear, &jmon, &jday, &jut);
3043 sgj = get_gregjul(gregflag, jyear);
3044 /*if ((eclflag = swe_lun_eclipse_how(t_ut, whicheph, geopos, attr, serr)) == ERR) {
3045 do_printf(serr);
3046 return ERR;
3047 }*/
3048 dt = (tret[3] - tret[2]) * 24 * 60;
3049 sprintf(s1, "%d min %4.2f sec", (int) dt, fmod(dt, 1) * 60);
3050 /* short output:
3051 * date, time of day, umbral magnitude, umbral duration, saros series, member number */
3052 sprintf(saros, "%d/%d", (int) attr[9], (int) attr[10]);
3053 sprintf(sout_short, "%s\t%2d.%2d.%4d%s\t%s\t%.3f\t%s\t%s\n", sout, jday, jmon, jyear, sgj, hms(jut,0), attr[8],s1, saros);
3054 sprintf(sout + strlen(sout), "%2d.%02d.%04d%s\t%s\t%.4f/%.4f\tsaros %s\t%.6f\n", jday, jmon, jyear, sgj, hms(jut,BIT_LZEROES), attr[0],attr[1], saros, t_ut);
3055 /* second line:
3056 * eclipse times, penumbral, partial, total begin and end */
3057 if (have_gap_parameter) sprintf(sout + strlen(sout), "\t");
3058 if (eclflag & SE_ECL_PENUMBBEG_VISIBLE)
3059 sprintf(sout + strlen(sout), " %s ", hms_from_tjd(tret[6]));
3060 else
3061 strcat(sout, " - ");
3062 if (have_gap_parameter) sprintf(sout + strlen(sout), "\t");
3063 if (eclflag & SE_ECL_PARTBEG_VISIBLE)
3064 sprintf(sout + strlen(sout), "%s ", hms_from_tjd(tret[2]));
3065 else
3066 strcat(sout, " - ");
3067 if (have_gap_parameter) sprintf(sout + strlen(sout), "\t");
3068 if (eclflag & SE_ECL_TOTBEG_VISIBLE)
3069 sprintf(sout + strlen(sout), "%s ", hms_from_tjd(tret[4]));
3070 else
3071 strcat(sout, " - ");
3072 if (have_gap_parameter) sprintf(sout + strlen(sout), "\t");
3073 if (eclflag & SE_ECL_TOTEND_VISIBLE)
3074 sprintf(sout + strlen(sout), "%s ", hms_from_tjd(tret[5]));
3075 else
3076 strcat(sout, " - ");
3077 if (have_gap_parameter) sprintf(sout + strlen(sout), "\t");
3078 if (eclflag & SE_ECL_PARTEND_VISIBLE)
3079 sprintf(sout + strlen(sout), "%s ", hms_from_tjd(tret[3]));
3080 else
3081 strcat(sout, " - ");
3082 if (have_gap_parameter) sprintf(sout + strlen(sout), "\t");
3083 if (eclflag & SE_ECL_PENUMBEND_VISIBLE)
3084 sprintf(sout + strlen(sout), "%s ", hms_from_tjd(tret[7]));
3085 else
3086 strcat(sout, " - ");
3087 if (have_gap_parameter) sprintf(sout + strlen(sout), "\t");
3088 sprintf(sout + strlen(sout), "dt=%.1f", swe_deltat_ex(tret[0], whicheph, serr) * 86400.0);
3089 strcat(sout, "\n");
3090 /* global lunar eclipse */
3091 } else {
3092 if ((eclflag = swe_lun_eclipse_when(t_ut, whicheph, search_flag,
3093 tret, direction_flag, serr)) == ERR) {
3094 do_printf(serr);
3095 return ERR;
3096 }
3097 t_ut = tret[0];
3098 if ((eclflag & SE_ECL_TOTAL)) {
3099 styp = "Total";
3100 strcpy(sout, "total ");
3101 ecl_type = ECL_LUN_TOTAL;
3102 }
3103 if ((eclflag & SE_ECL_PENUMBRAL)) {
3104 styp = "Penumbral";
3105 strcpy(sout, "penumb. ");
3106 ecl_type = ECL_LUN_PENUMBRAL;
3107 }
3108 if ((eclflag & SE_ECL_PARTIAL)) {
3109 styp = "Partial";
3110 strcpy(sout, "partial ");
3111 ecl_type = ECL_LUN_PARTIAL;
3112 }
3113 strcat(sout, "lunar eclipse\t");
3114 if ((eclflag = swe_lun_eclipse_how(t_ut, whicheph, geopos, attr, serr)) ==
3115 ERR) {
3116 do_printf(serr);
3117 return ERR;
3118 }
3119 if (time_flag & (BIT_TIME_LMT | BIT_TIME_LAT)) {
3120 for (i = 0; i < 10; i++) {
3121 if (tret[i] != 0) {
3122 retc = ut_to_lmt_lat(tret[i], geopos, &(tret[i]), serr);
3123 if (retc == ERR) {
3124 do_printf(serr);
3125 return ERR;
3126 }
3127 }
3128 }
3129 }
3130 t_ut = tret[0];
3131 rval = swe_calc_ut(t_ut, SE_MOON, whicheph|SEFLG_EQUATORIAL, xx, s1);
3132 if (rval < 0)
3133 strcat(s1, "\n");
3134 do_printf(s1);
3135 swe_revjul(t_ut, gregflag, &jyear, &jmon, &jday, &jut);
3136 geopos_max[0] = swe_degnorm(xx[0]- swe_sidtime(t_ut) * 15);
3137 if (geopos_max[0] > 180) geopos_max[0] -= 360;
3138 geopos_max[1] = xx[1];
3139 sgj = get_gregjul(gregflag, jyear);
3140 dt = (tret[3] - tret[2]) * 24 * 60;
3141 sprintf(s1, "%d min %4.2f sec", (int) dt, fmod(dt, 1) * 60);
3142 /* short output:
3143 * date, time of day, umbral magnitude, umbral duration, saros series, member number */
3144 sprintf(saros, "%d/%d", (int) attr[9], (int) attr[10]);
3145 sprintf(sout_short, "%s\t%2d.%2d.%4d%s\t%s\t%.3f\t%s\t%s\n", sout, jday, jmon, jyear, sgj, hms(jut,0), attr[8],s1, saros);
3146 //sprintf(sout + strlen(sout), "%2d.%02d.%04d%s\t%s\t%.4f/%.4f\tsaros %s\t%.6f\tdt=%.2f\n", jday, jmon, jyear, sgj, hms(jut,BIT_LZEROES), attr[0],attr[1], saros, t_ut, swe_deltat_ex(t_ut, whicheph, serr) * 86400);
3147 sprintf(sout + strlen(sout), "%2d.%02d.%04d%s\t%s\t%.4f/%.4f\tsaros %s\t%.6f\n", jday, jmon, jyear, sgj, hms(jut,BIT_LZEROES), attr[0],attr[1], saros, t_ut);
3148 /* second line:
3149 * eclipse times, penumbral, partial, total begin and end */
3150 if (have_gap_parameter) sprintf(sout + strlen(sout), "\t");
3151 sprintf(sout + strlen(sout), " %s ", hms_from_tjd(tret[6]));
3152 if (have_gap_parameter) sprintf(sout + strlen(sout), "\t");
3153 if (tret[2] != 0)
3154 sprintf(sout + strlen(sout), "%s ", hms_from_tjd(tret[2]));
3155 else
3156 strcat(sout, " - ");
3157 if (have_gap_parameter) sprintf(sout + strlen(sout), "\t");
3158 if (tret[4] != 0)
3159 sprintf(sout + strlen(sout), "%s ", hms_from_tjd(tret[4]));
3160 else
3161 strcat(sout, " - ");
3162 if (have_gap_parameter) sprintf(sout + strlen(sout), "\t");
3163 if (tret[5] != 0)
3164 sprintf(sout + strlen(sout), "%s ", hms_from_tjd(tret[5]));
3165 else
3166 strcat(sout, " - ");
3167 if (have_gap_parameter) sprintf(sout + strlen(sout), "\t");
3168 if (tret[3] != 0)
3169 sprintf(sout + strlen(sout), "%s ", hms_from_tjd(tret[3]));
3170 else
3171 strcat(sout, " - ");
3172 if (have_gap_parameter) sprintf(sout + strlen(sout), "\t");
3173 sprintf(sout + strlen(sout), "%s", hms_from_tjd(tret[7]));
3174 if (have_gap_parameter) sprintf(sout + strlen(sout), "\t");
3175 sprintf(sout + strlen(sout), "dt=%.1f", swe_deltat_ex(tret[0], whicheph, serr) * 86400.0);
3176 strcat(sout, "\n");
3177 if (special_mode & SP_MODE_HOCAL) {
3178 swe_split_deg(jut, SE_SPLIT_DEG_ROUND_MIN, &ihou, &imin, &isec, &dfrc, &isgn);
3179 sprintf(sout, "\"%04d%s %02d %02d %02d.%02d %d\",\n", jyear, sgj, jmon, jday, ihou, imin, ecl_type);
3180 }
3181 sprintf(sout + strlen(sout), "\t%s\t%s\n", strcpy(s1, dms(geopos_max[0], BIT_ROUND_SEC)), strcpy(s2, dms(geopos_max[1], BIT_ROUND_SEC)));
3182 }
3183 //dt = (tret[7] - tret[6]) * 24 * 60;
3184 //sprintf(sout + strlen(sout), "\t%d min %4.2f sec\n", (int) dt, fmod(dt, 1) * 60);
3185 if (have_gap_parameter) insert_gap_string_for_tabs(sout, gap);
3186 if (short_output) {
3187 do_printf(sout_short);
3188 } else {
3189 do_printf(sout);
3190 }
3191 if (with_chart_link) {
3192 char snat[AS_MAXCH];
3193 char stim[AS_MAXCH];
3194 int iflg = 0;
3195 char cal = gregflag ? 'g' : 'j';
3196 strcpy(stim, hms(jut,BIT_LZEROES));
3197 format_lon_lat(slon, slat, geopos_max[0], geopos_max[1]);
3198 while (*stim == ' ') our_strcpy(stim, stim + 1);
3199 if (*stim == '0') our_strcpy(stim, stim + 1);
3200 sprintf(snat, "Lunar Eclipse %s,%s,e,%d,%d,%d,%s,h0e,%cnu,%d,Moon Zenith location,,%s,%s,u,0,0,0", saros, styp, jday, jmon, jyear, stim, cal, iflg, slon, slat);
3201 sprintf(sout, "<a href='https://www.astro.com/cgi/chart.cgi?muasp=1;nhor=1;act=chmnat;nd1=%s;rs=1;iseclipse=1' target='eclipse'>chart link</a>\n\n", snat);
3202 do_printf(sout);
3203 }
3204 }
3205 if (with_chart_link) do_printf("</pre>\n");
3206 return OK;
3207 }
3208
call_solar_eclipse(double t_ut,int32 whicheph,int32 special_mode,double * geopos,char * serr)3209 static int32 call_solar_eclipse(double t_ut, int32 whicheph, int32 special_mode, double *geopos, char *serr)
3210 {
3211 int i, ii, retc = OK, eclflag, ecl_type = 0;
3212 double dt, tret[30], attr[30], geopos_max[3];
3213 char slon[8], slat[8], saros[20];
3214 char s1[AS_MAXCH], s2[AS_MAXCH], sout_short[AS_MAXCH], *styp = "none", *sgj;
3215 AS_BOOL has_found = FALSE;
3216 /* no selective eclipse type set, set all */
3217 if (with_chart_link) do_printf("<pre>");
3218 if ((search_flag & SE_ECL_ALLTYPES_SOLAR) == 0)
3219 search_flag |= SE_ECL_ALLTYPES_SOLAR;
3220 /* for local eclipses: set geographic position of observer */
3221 if (special_mode & SP_MODE_LOCAL) {
3222 swe_set_topo(geopos[0], geopos[1], geopos[2]);
3223 // "geo. long 8.000000, lat 47.000000, alt 0.000000"
3224 if (with_header)
3225 printf("\ngeo. long %f, lat %f, alt %f", geopos[0], geopos[1], geopos[2]);
3226 }
3227 do_printf("\n");
3228 for (ii = 0; ii < nstep; ii++, t_ut += direction) {
3229 *sout = '\0';
3230 /* swetest -solecl -local -geopos...
3231 * find next solar eclipse observable from a given geographic position */
3232 if (special_mode & SP_MODE_LOCAL) {
3233 if ((eclflag = swe_sol_eclipse_when_loc(t_ut, whicheph, geopos, tret,
3234 attr, direction_flag, serr)) == ERR) {
3235 do_printf(serr);
3236 return ERR;
3237 } else {
3238 has_found = FALSE;
3239 t_ut = tret[0];
3240 if ((search_flag & SE_ECL_TOTAL) && (eclflag & SE_ECL_TOTAL)) {
3241 strcpy(sout, "total ");
3242 has_found = TRUE;
3243 ecl_type = ECL_SOL_TOTAL;
3244 }
3245 if ((search_flag & SE_ECL_ANNULAR) && (eclflag & SE_ECL_ANNULAR)) {
3246 strcpy(sout, "annular ");
3247 has_found = TRUE;
3248 ecl_type = ECL_SOL_ANNULAR;
3249 }
3250 if ((search_flag & SE_ECL_PARTIAL) && (eclflag & SE_ECL_PARTIAL)) {
3251 strcpy(sout, "partial ");
3252 has_found = TRUE;
3253 ecl_type = ECL_SOL_PARTIAL;
3254 }
3255 if (have_gap_parameter) sprintf(sout + strlen(sout), "\t");
3256 if (!has_found) {
3257 ii--;
3258 } else {
3259 swe_calc(t_ut + swe_deltat_ex(t_ut, whicheph, serr), SE_ECL_NUT, 0, x, serr);
3260 if (time_flag & (BIT_TIME_LMT | BIT_TIME_LAT)) {
3261 for (i = 0; i < 10; i++) {
3262 if (tret[i] != 0) {
3263 retc = ut_to_lmt_lat(tret[i], geopos, &(tret[i]), serr);
3264 if (retc == ERR) {
3265 do_printf(serr);
3266 return ERR;
3267 }
3268 }
3269 }
3270 }
3271 t_ut = tret[0];
3272 swe_revjul(t_ut, gregflag, &jyear, &jmon, &jday, &jut);
3273 dt = (tret[3] - tret[2]) * 24 * 60;
3274 sgj = get_gregjul(gregflag, jyear);
3275 sprintf(saros, "%d/%d", (int) attr[9], (int) attr[10]);
3276 sprintf(sout + strlen(sout), "%2d.%02d.%04d%s\t%s\t%.4f/%.4f/%.4f\tsaros %s\t%.6f\n", jday, jmon, jyear, sgj,hms(jut,BIT_LZEROES), attr[8], attr[0], attr[2], saros, t_ut);
3277 sprintf(sout + strlen(sout), "\t%d min %4.2f sec\t", (int) dt, fmod(dt, 1) * 60);
3278 if (eclflag & SE_ECL_1ST_VISIBLE) {
3279 sprintf(sout + strlen(sout), "%s ", hms_from_tjd(tret[1]));
3280 } else {
3281 strcat(sout, " - ");
3282 }
3283 if (have_gap_parameter) sprintf(sout + strlen(sout), "\t");
3284 if (eclflag & SE_ECL_2ND_VISIBLE) {
3285 sprintf(sout + strlen(sout), "%s ", hms_from_tjd(tret[2]));
3286 } else {
3287 strcat(sout, " - ");
3288 }
3289 if (have_gap_parameter) sprintf(sout + strlen(sout), "\t");
3290 if (eclflag & SE_ECL_3RD_VISIBLE) {
3291 sprintf(sout + strlen(sout), "%s ", hms_from_tjd(tret[3]));
3292 } else {
3293 strcat(sout, " - ");
3294 }
3295 if (have_gap_parameter) sprintf(sout + strlen(sout), "\t");
3296 if (eclflag & SE_ECL_4TH_VISIBLE) {
3297 sprintf(sout + strlen(sout), "%s ", hms_from_tjd(tret[4]));
3298 } else {
3299 strcat(sout, " - ");
3300 }
3301 if (have_gap_parameter) sprintf(sout + strlen(sout), "\t");
3302 #if 0
3303 sprintf(sout + strlen(sout), "\t%d min %4.2f sec %s %s %s %s",
3304 (int) dt, fmod(dt, 1) * 60,
3305 strcpy(s1, hms(fmod(tret[1] + 0.5, 1) * 24, BIT_LZEROES)),
3306 strcpy(s3, hms(fmod(tret[2] + 0.5, 1) * 24, BIT_LZEROES)),
3307 strcpy(s4, hms(fmod(tret[3] + 0.5, 1) * 24, BIT_LZEROES)),
3308 strcpy(s2, hms(fmod(tret[4] + 0.5, 1) * 24, BIT_LZEROES)));
3309 #endif
3310 sprintf(sout + strlen(sout), "dt=%.1f", swe_deltat_ex(tret[0], whicheph, serr) * 86400.0);
3311 strcat(sout, "\n");
3312 if (have_gap_parameter) insert_gap_string_for_tabs(sout, gap);
3313 do_printf(sout);
3314 }
3315 }
3316 } /* endif search_local */
3317 /* swetest -solecl
3318 * find next solar eclipse observable from anywhere on earth */
3319 if (!(special_mode & SP_MODE_LOCAL)) {
3320 if ((eclflag = swe_sol_eclipse_when_glob(t_ut, whicheph, search_flag,
3321 tret, direction_flag, serr)) == ERR) {
3322 do_printf(serr);
3323 return ERR;
3324 }
3325 t_ut = tret[0];
3326 if ((eclflag & SE_ECL_TOTAL)) {
3327 styp = "Total";
3328 strcpy(sout, "total");
3329 ecl_type = ECL_SOL_TOTAL;
3330 }
3331 if ((eclflag & SE_ECL_ANNULAR)) {
3332 styp = "Annular";
3333 strcpy(sout, "annular");
3334 ecl_type = ECL_SOL_ANNULAR;
3335 }
3336 if ((eclflag & SE_ECL_ANNULAR_TOTAL)) {
3337 styp = "Annular-Total";
3338 strcpy(sout, "ann-tot");
3339 ecl_type = ECL_SOL_ANNULAR; /* by Alois: what is this ? */
3340 }
3341 if ((eclflag & SE_ECL_PARTIAL)) {
3342 styp = "Partial";
3343 strcpy(sout, "partial");
3344 ecl_type = ECL_SOL_PARTIAL;
3345 }
3346 if ((eclflag & SE_ECL_NONCENTRAL) && !(eclflag & SE_ECL_PARTIAL))
3347 strcat(sout, " non-central");
3348 sprintf(sout + strlen(sout), " solar\t");
3349 swe_sol_eclipse_where(t_ut, whicheph, geopos_max, attr, serr);
3350 if (time_flag & (BIT_TIME_LMT | BIT_TIME_LAT)) {
3351 for (i = 0; i < 10; i++) {
3352 if (tret[i] != 0) {
3353 retc = ut_to_lmt_lat(tret[i], geopos, &(tret[i]), serr);
3354 if (retc == ERR) {
3355 do_printf(serr);
3356 return ERR;
3357 }
3358 }
3359 }
3360 }
3361 swe_revjul(tret[0], gregflag, &jyear, &jmon, &jday, &jut);
3362 sgj = get_gregjul(gregflag, jyear);
3363 sprintf(saros, "%d/%d", (int) attr[9], (int) attr[10]);
3364 sprintf(sout_short, "%s\t%2d.%2d.%4d%s\t%s\t%.3f", sout, jday, jmon, jyear, sgj, hms(jut,0), attr[8]);
3365 sprintf(sout + strlen(sout), "%2d.%02d.%04d%s\t%s\t%f km\t%.4f/%.4f/%.4f\tsaros %s\t%.6f\n", jday, jmon, jyear, sgj, hms(jut,0), attr[3], attr[8], attr[0], attr[2], saros, tret[0]);
3366 sprintf(sout + strlen(sout), "\t%s ", hms_from_tjd(tret[2]));
3367 if (have_gap_parameter) sprintf(sout + strlen(sout), "\t");
3368 if (tret[4] != 0) {
3369 sprintf(sout + strlen(sout), "%s ", hms_from_tjd(tret[4]));
3370 } else {
3371 strcat(sout, " - ");
3372 }
3373 if (have_gap_parameter) sprintf(sout + strlen(sout), "\t");
3374 if (tret[5] != 0) {
3375 sprintf(sout + strlen(sout), "%s ", hms_from_tjd(tret[5]));
3376 } else {
3377 strcat(sout, " - ");
3378 }
3379 if (have_gap_parameter) sprintf(sout + strlen(sout), "\t");
3380 sprintf(sout + strlen(sout), "%s", hms_from_tjd(tret[3]));
3381 if (have_gap_parameter) sprintf(sout + strlen(sout), "\t");
3382 sprintf(sout + strlen(sout), "dt=%.1f", swe_deltat_ex(tret[0], whicheph, serr) * 86400.0);
3383 strcat(sout, "\n");
3384 sprintf(sout + strlen(sout), "\t%s\t%s", strcpy(s1, dms(geopos_max[0], BIT_ROUND_SEC)), strcpy(s2, dms(geopos_max[1], BIT_ROUND_SEC)));
3385 strcat(sout, "\t");
3386 strcat(sout_short, "\t");
3387 if (!(eclflag & SE_ECL_PARTIAL) && !(eclflag & SE_ECL_NONCENTRAL)) {
3388 if ((eclflag = swe_sol_eclipse_when_loc(t_ut - 10, whicheph, geopos_max, tret, attr, 0, serr)) == ERR) {
3389 do_printf(serr);
3390 return ERR;
3391 }
3392 if (fabs(tret[0] - t_ut) > 2) {
3393 do_printf("when_loc returns wrong date\n");
3394 }
3395 dt = (tret[3] - tret[2]) * 24 * 60;
3396 sprintf(s1, "%d min %4.2f sec", (int) dt, fmod(dt, 1) * 60);
3397 strcat(sout, s1);
3398 strcat(sout_short, s1);
3399 }
3400 sprintf(sout_short + strlen(sout_short), "\t%d\t%d", (int) attr[9], (int) attr[10]);
3401 strcat(sout, "\n");
3402 strcat(sout_short, "\n");
3403 if (special_mode & SP_MODE_HOCAL) {
3404 int ihou, imin, isec, isgn;
3405 double dfrc;
3406 swe_split_deg(jut, SE_SPLIT_DEG_ROUND_MIN, &ihou, &imin, &isec, &dfrc, &isgn);
3407 sprintf(sout, "\"%04d%s %02d %02d %02d.%02d %d\",\n", jyear, sgj, jmon, jday, ihou, imin, ecl_type);
3408 }
3409 /*printf("len=%ld\n", strlen(sout));*/
3410 if (short_output) {
3411 do_printf(sout_short);
3412 } else {
3413 if (have_gap_parameter) insert_gap_string_for_tabs(sout, gap);
3414 do_printf(sout);
3415 }
3416 if (with_chart_link) {
3417 char snat[AS_MAXCH];
3418 char stim[AS_MAXCH];
3419 int iflg = 0; // NAT_IFLG_UNKNOWN_TIME;
3420 char cal = gregflag ? 'g' : 'j';
3421 format_lon_lat(slon, slat, geopos_max[0], geopos_max[1]);
3422 strcpy(stim, hms(jut,BIT_LZEROES));
3423 while (*stim == ' ') our_strcpy(stim, stim + 1);
3424 if (*stim == '0') our_strcpy(stim, stim + 1);
3425 sprintf(snat, "Solar Eclipse %s,%s,e,%d,%d,%d,%s,h0e,%cnu,%d,Location of Maximum,,%s,%s,u,0,0,0", saros, styp, jday, jmon, jyear, stim, cal, iflg, slon, slat);
3426 sprintf(sout, "<a href='https://www.astro.com/cgi/chart.cgi?muasp=1;nhor=1;act=chmnat;nd1=%s;rs=1;iseclipse=1;topo=1' target='eclipse'>chart link</a>\n\n", snat);
3427 do_printf(sout);
3428 }
3429 }
3430 }
3431 if (with_chart_link) do_printf("</pre>\n");
3432 return OK;
3433 }
3434
call_lunar_occultation(double t_ut,int32 ipl,char * star,int32 whicheph,int32 special_mode,double * geopos,char * serr)3435 static int32 call_lunar_occultation(double t_ut, int32 ipl, char *star, int32 whicheph, int32 special_mode, double *geopos, char *serr)
3436 {
3437 int i, ii, ecl_type = 0, eclflag, retc = OK;
3438 double dt, tret[30], attr[30], geopos_max[3];
3439 char s1[AS_MAXCH], s2[AS_MAXCH];
3440 AS_BOOL has_found = FALSE;
3441 int nloops = 0;
3442 /* no selective eclipse type set, set all */
3443 if ((search_flag & SE_ECL_ALLTYPES_SOLAR) == 0)
3444 search_flag |= SE_ECL_ALLTYPES_SOLAR;
3445 /* for local occultations: set geographic position of observer */
3446 if (special_mode & SP_MODE_LOCAL) {
3447 swe_set_topo(geopos[0], geopos[1], geopos[2]);
3448 if (with_header)
3449 printf("\ngeo. long %f, lat %f, alt %f", geopos[0], geopos[1], geopos[2]);
3450 }
3451 do_printf("\n");
3452 for (ii = 0; ii < nstep; ii++) {
3453 *sout = '\0';
3454 nloops++;
3455 if (nloops > SEARCH_RANGE_LUNAR_CYCLES) {
3456 sprintf(serr, "event search ended after %d lunar cycles at jd=%f\n", SEARCH_RANGE_LUNAR_CYCLES, t_ut);
3457 do_printf(serr);
3458 return ERR;
3459 }
3460 if (special_mode & SP_MODE_LOCAL) {
3461 /* * local search for occultation, test one lunar cycle only (SE_ECL_ONE_TRY) */
3462 if (ipl != SE_SUN) {
3463 search_flag &= ~(SE_ECL_ANNULAR|SE_ECL_ANNULAR_TOTAL);
3464 if (search_flag == 0)
3465 search_flag = SE_ECL_ALLTYPES_SOLAR;
3466 }
3467 if ((eclflag = swe_lun_occult_when_loc(t_ut, ipl, star, whicheph, geopos, tret, attr, direction_flag|SE_ECL_ONE_TRY, serr)) == ERR) {
3468 do_printf(serr);
3469 return ERR;
3470 } else if (eclflag == 0) { /* event not found, try next conjunction */
3471 t_ut = tret[0] + direction * 10; /* try again with start date increased by 10 */
3472 ii--;
3473 } else {
3474 t_ut = tret[0];
3475 if (time_flag & (BIT_TIME_LMT | BIT_TIME_LAT)) {
3476 for (i = 0; i < 10; i++) {
3477 if (tret[i] != 0) {
3478 retc = ut_to_lmt_lat(tret[i], geopos, &(tret[i]), serr);
3479 if (retc == ERR) {
3480 do_printf(serr);
3481 return ERR;
3482 }
3483 }
3484 }
3485 }
3486 has_found = FALSE;
3487 *sout = '\0';
3488 if ((search_flag & SE_ECL_TOTAL) && (eclflag & SE_ECL_TOTAL)) {
3489 strcat(sout, "total");
3490 has_found = TRUE;
3491 ecl_type = ECL_SOL_TOTAL;
3492 }
3493 if ((search_flag & SE_ECL_ANNULAR) && (eclflag & SE_ECL_ANNULAR)) {
3494 strcat(sout, "annular");
3495 has_found = TRUE;
3496 ecl_type = ECL_SOL_ANNULAR;
3497 }
3498 if ((search_flag & SE_ECL_PARTIAL) && (eclflag & SE_ECL_PARTIAL)) {
3499 strcat(sout, "partial");
3500 has_found = TRUE;
3501 ecl_type = ECL_SOL_PARTIAL;
3502 }
3503 if (ipl != SE_SUN) {
3504 if ((eclflag & SE_ECL_OCC_BEG_DAYLIGHT) && (eclflag & SE_ECL_OCC_END_DAYLIGHT))
3505 strcat(sout, "(daytime)"); /* occultation occurs during the day */
3506 else if (eclflag & SE_ECL_OCC_BEG_DAYLIGHT)
3507 strcat(sout, "(sunset) "); /* occultation occurs during the day */
3508 else if (eclflag & SE_ECL_OCC_END_DAYLIGHT)
3509 strcat(sout, "(sunrise)"); /* occultation occurs during the day */
3510 }
3511 if (have_gap_parameter) sprintf(sout + strlen(sout), "\t");
3512 while (strlen(sout) < 17)
3513 strcat(sout, " ");
3514 if (!has_found) {
3515 ii--;
3516 } else {
3517 swe_calc_ut(t_ut, SE_ECL_NUT, 0, x, serr);
3518 swe_revjul(tret[0], gregflag, &jyear, &jmon, &jday, &jut);
3519 dt = (tret[3] - tret[2]) * 24 * 60;
3520 sprintf(sout + strlen(sout), "%2d.%02d.%04d\t%s\t%f\t%.6f\n", jday, jmon, jyear, hms(jut,BIT_LZEROES), attr[0], tret[0]);
3521 sprintf(sout + strlen(sout), "\t%d min %4.2f sec\t", (int) dt, fmod(dt, 1) * 60);
3522 if (eclflag & SE_ECL_1ST_VISIBLE)
3523 sprintf(sout + strlen(sout), "%s ", hms_from_tjd(tret[1]));
3524 else
3525 strcat(sout, " - ");
3526 if (have_gap_parameter) sprintf(sout + strlen(sout), "\t");
3527 if (eclflag & SE_ECL_2ND_VISIBLE)
3528 sprintf(sout + strlen(sout), "%s ", hms_from_tjd(tret[2]));
3529 else
3530 strcat(sout, " - ");
3531 if (have_gap_parameter) sprintf(sout + strlen(sout), "\t");
3532 if (eclflag & SE_ECL_3RD_VISIBLE)
3533 sprintf(sout + strlen(sout), "%s ", hms_from_tjd(tret[3]));
3534 else
3535 strcat(sout, " - ");
3536 if (have_gap_parameter) sprintf(sout + strlen(sout), "\t");
3537 if (eclflag & SE_ECL_4TH_VISIBLE)
3538 sprintf(sout + strlen(sout), "%s ", hms_from_tjd(tret[4]));
3539 else
3540 strcat(sout, " - ");
3541 if (have_gap_parameter) sprintf(sout + strlen(sout), "\t");
3542 #if 0
3543 sprintf(sout + strlen(sout), "\t%d min %4.2f sec %s %s %s %s",
3544 (int) dt, fmod(dt, 1) * 60,
3545 strcpy(s1, hms(fmod(tret[1] + 0.5, 1) * 24, BIT_LZEROES)),
3546 strcpy(s3, hms(fmod(tret[2] + 0.5, 1) * 24, BIT_LZEROES)),
3547 strcpy(s4, hms(fmod(tret[3] + 0.5, 1) * 24, BIT_LZEROES)),
3548 strcpy(s2, hms(fmod(tret[4] + 0.5, 1) * 24, BIT_LZEROES)));
3549 #endif
3550 sprintf(sout + strlen(sout), "dt=%.1f", swe_deltat_ex(tret[0], whicheph, serr) * 86400.0);
3551 strcat(sout, "\n");
3552 if (have_gap_parameter) insert_gap_string_for_tabs(sout, gap);
3553 do_printf(sout);
3554 }
3555 }
3556 } /* endif search_local */
3557 if (!(special_mode & SP_MODE_LOCAL)) {
3558 /* * global search for occultations, test one lunar cycle only (SE_ECL_ONE_TRY) */
3559 if ((eclflag = swe_lun_occult_when_glob(t_ut, ipl, star, whicheph, search_flag, tret, direction_flag|SE_ECL_ONE_TRY, serr)) == ERR) {
3560 do_printf(serr);
3561 return ERR;
3562 }
3563 if (eclflag == 0) { /* no occltation was found at next conjunction, try next conjunction */
3564 t_ut = tret[0] + direction;
3565 ii--;
3566 continue;
3567 }
3568 if ((eclflag & SE_ECL_TOTAL)) {
3569 strcpy(sout, "total ");
3570 ecl_type = ECL_SOL_TOTAL;
3571 }
3572 if ((eclflag & SE_ECL_ANNULAR)) {
3573 strcpy(sout, "annular ");
3574 ecl_type = ECL_SOL_ANNULAR;
3575 }
3576 if ((eclflag & SE_ECL_ANNULAR_TOTAL)) {
3577 strcpy(sout, "ann-tot ");
3578 ecl_type = ECL_SOL_ANNULAR; /* by Alois: what is this ? */
3579 }
3580 if ((eclflag & SE_ECL_PARTIAL)) {
3581 strcpy(sout, "partial ");
3582 ecl_type = ECL_SOL_PARTIAL;
3583 }
3584 if ((eclflag & SE_ECL_NONCENTRAL) && !(eclflag & SE_ECL_PARTIAL))
3585 strcat(sout, "non-central ");
3586 t_ut = tret[0];
3587 swe_lun_occult_where(t_ut, ipl, star, whicheph, geopos_max, attr, serr);
3588 /* for (i = 0; i < 8; i++) {
3589 printf("attr[%d]=%.17f\n", i, attr[i]);
3590 } */
3591 if (time_flag & (BIT_TIME_LMT | BIT_TIME_LAT)) {
3592 for (i = 0; i < 10; i++) {
3593 if (tret[i] != 0) {
3594 retc = ut_to_lmt_lat(tret[i], geopos, &(tret[i]), serr);
3595 if (retc == ERR) {
3596 do_printf(serr);
3597 return ERR;
3598 }
3599 }
3600 }
3601 }
3602 swe_revjul(tret[0], gregflag, &jyear, &jmon, &jday, &jut);
3603 sprintf(sout + strlen(sout), "%2d.%02d.%04d\t%s\t%f km\t%f\t%.6f\n", jday, jmon, jyear, hms(jut,BIT_LZEROES), attr[3], attr[0], tret[0]);
3604 sprintf(sout + strlen(sout), "\t%s ", hms_from_tjd(tret[2]));
3605 if (have_gap_parameter) sprintf(sout + strlen(sout), "\t");
3606 if (tret[4] != 0)
3607 sprintf(sout + strlen(sout), "%s ", hms_from_tjd(tret[4]));
3608 else
3609 strcat(sout, " - ");
3610 if (have_gap_parameter) sprintf(sout + strlen(sout), "\t");
3611 if (tret[5] != 0)
3612 sprintf(sout + strlen(sout), "%s ", hms_from_tjd(tret[5]));
3613 else
3614 strcat(sout, " - ");
3615 if (have_gap_parameter) sprintf(sout + strlen(sout), "\t");
3616 sprintf(sout + strlen(sout), "%s", hms_from_tjd(tret[3]));
3617 if (have_gap_parameter) sprintf(sout + strlen(sout), "\t");
3618 sprintf(sout + strlen(sout), "dt=%.1f", swe_deltat_ex(tret[0], whicheph, serr) * 86400.0);
3619 strcat(sout, "\n");
3620 sprintf(sout + strlen(sout), "\t%s\t%s", strcpy(s1, dms(geopos_max[0], BIT_ROUND_MIN)), strcpy(s2, dms(geopos_max[1], BIT_ROUND_MIN)));
3621 if (!(eclflag & SE_ECL_PARTIAL) && !(eclflag & SE_ECL_NONCENTRAL)) {
3622 if ((eclflag = swe_lun_occult_when_loc(t_ut - 10, ipl, star, whicheph, geopos_max, tret, attr, 0, serr)) == ERR) {
3623 do_printf(serr);
3624 return ERR;
3625 }
3626 if (fabs(tret[0] - t_ut) > 2)
3627 do_printf("when_loc returns wrong date\n");
3628 dt = (tret[3] - tret[2]) * 24 * 60;
3629 sprintf(sout + strlen(sout), "\t%d min %4.2f sec", (int) dt, fmod(dt, 1) * 60);
3630 }
3631 strcat(sout, "\n");
3632 if (have_gap_parameter) insert_gap_string_for_tabs(sout, gap);
3633 if (special_mode & SP_MODE_HOCAL) {
3634 int ihou, imin, isec, isgn;
3635 double dfrc;
3636 swe_split_deg(jut, SE_SPLIT_DEG_ROUND_MIN, &ihou, &imin, &isec, &dfrc, &isgn);
3637 sprintf(sout, "\"%04d %02d %02d %02d.%02d %d\",\n", jyear, jmon, jday, ihou, imin, ecl_type);
3638 }
3639 do_printf(sout);
3640 }
3641 t_ut += direction;
3642 }
3643 return OK;
3644 }
3645
do_print_heliacal(double * dret,int32 event_type,char * obj_name)3646 static void do_print_heliacal(double *dret, int32 event_type, char *obj_name)
3647 {
3648 char *sevtname[] = {"",
3649 "heliacal rising ",
3650 "heliacal setting",
3651 "evening first ",
3652 "morning last ",
3653 "evening rising ",
3654 "morning setting ",};
3655 char *stz = "UT";
3656 char stim0[40], stim1[40], stim2[40];
3657 if (time_flag & BIT_TIME_LMT)
3658 stz = "LMT";
3659 if (time_flag & BIT_TIME_LAT)
3660 stz = "LAT";
3661 *sout = '\0';
3662 swe_revjul(dret[0], gregflag, &jyear, &jmon, &jday, &jut);
3663 if (event_type <= 4) {
3664 if (hel_using_AV) {
3665 strcpy(stim0, hms_from_tjd(dret[0]));
3666 remove_whitespace(stim0);
3667 /* The following line displays only the beginning of visibility. */
3668 sprintf(sout + strlen(sout), "%s %s: %d/%02d/%02d %s %s (%.5f)\n", obj_name, sevtname[event_type], jyear, jmon, jday, stim0, stz, dret[0]);
3669 } else {
3670 /* display the moment of beginning and optimum visibility */
3671 strcpy(stim0, hms_from_tjd(dret[0]));
3672 strcpy(stim1, hms_from_tjd(dret[1]));
3673 strcpy(stim2, hms_from_tjd(dret[2]));
3674 remove_whitespace(stim0);
3675 remove_whitespace(stim1);
3676 remove_whitespace(stim2);
3677 sprintf(sout + strlen(sout), "%s %s: %d/%02d/%02d %s %s (%.5f), opt %s, end %s, dur %.1f min\n", obj_name, sevtname[event_type], jyear, jmon, jday, stim0, stz, dret[0], stim1, stim2, (dret[2] - dret[0]) * 1440);
3678 }
3679 } else {
3680 strcpy(stim0, hms_from_tjd(dret[0]));
3681 remove_whitespace(stim0);
3682 sprintf(sout + strlen(sout), "%s %s: %d/%02d/%02d %s %s (%f)\n", obj_name, sevtname[event_type], jyear, jmon, jday, stim0, stz, dret[0]);
3683 }
3684 do_printf(sout);
3685 }
3686
call_heliacal_event(double t_ut,int32 ipl,char * star,int32 whicheph,double * geopos,double * datm,double * dobs,char * serr)3687 static int32 call_heliacal_event(double t_ut, int32 ipl, char *star, int32 whicheph, double *geopos, double *datm, double *dobs, char *serr)
3688 {
3689 int ii, retc = OK, event_type = 0, retflag;
3690 double dret[40], tsave1, tsave2 = 0;
3691 char obj_name[AS_MAXCH];
3692 helflag |= whicheph;
3693 /* if invalid heliacal event type was required, set 0 for any type */
3694 if (search_flag < 0 || search_flag > 6)
3695 search_flag = 0;
3696 /* optical instruments used: */
3697 if (dobs[3] > 0)
3698 helflag |= SE_HELFLAG_OPTICAL_PARAMS;
3699 if (hel_using_AV)
3700 helflag |= SE_HELFLAG_AV;
3701 if (ipl == SE_FIXSTAR)
3702 strcpy(obj_name, star);
3703 else
3704 swe_get_planet_name(ipl, obj_name);
3705 if (with_header) {
3706 printf("\ngeo. long %f, lat %f, alt %f", geopos[0], geopos[1], geopos[2]);
3707 do_printf("\n");
3708 }
3709 for (ii = 0; ii < nstep; ii++, t_ut = dret[0] + 1) {
3710 *sout = '\0';
3711 if (search_flag > 0)
3712 event_type = search_flag;
3713 else if (ipl == SE_MOON)
3714 event_type = SE_EVENING_FIRST;
3715 else
3716 event_type = SE_HELIACAL_RISING;
3717 retflag = swe_heliacal_ut(t_ut, geopos, datm, dobs, obj_name, event_type, helflag, dret, serr);
3718 if (retflag == ERR) {
3719 do_printf(serr);
3720 return ERR;
3721 }
3722 if (time_flag & (BIT_TIME_LMT | BIT_TIME_LAT)) {
3723 retc = ut_to_lmt_lat(dret[0], geopos, &(dret[0]), serr);
3724 if (retc != ERR) retc = ut_to_lmt_lat(dret[1], geopos, &(dret[1]), serr);
3725 if (retc != ERR) retc = ut_to_lmt_lat(dret[2], geopos, &(dret[2]), serr);
3726 if (retc == ERR) {
3727 do_printf(serr);
3728 return ERR;
3729 }
3730 }
3731 do_print_heliacal(dret, event_type, obj_name);
3732 /* list all events within synodic cycle */
3733 if (search_flag == 0) {
3734 if (ipl == SE_VENUS || ipl == SE_MERCURY) {
3735 /* we have heliacal rising (morning first), now find morning last */
3736 event_type = SE_MORNING_LAST;
3737 retflag = swe_heliacal_ut(dret[0], geopos, datm, dobs, obj_name, event_type, helflag, dret, serr);
3738 if (retflag == ERR) {
3739 do_printf(serr);
3740 return ERR;
3741 }
3742 if (time_flag & (BIT_TIME_LMT | BIT_TIME_LAT)) {
3743 retc = ut_to_lmt_lat(dret[0], geopos, &(dret[0]), serr);
3744 if (retc != ERR) retc = ut_to_lmt_lat(dret[1], geopos, &(dret[1]), serr);
3745 if (retc != ERR) retc = ut_to_lmt_lat(dret[2], geopos, &(dret[2]), serr);
3746 if (retc == ERR) {
3747 do_printf(serr);
3748 return ERR;
3749 }
3750 }
3751 do_print_heliacal(dret, event_type, obj_name);
3752 tsave1 = dret[0];
3753 /* mercury can have several evening appearances without any morning
3754 * appearances in betweeen. We have to find out when the next
3755 * morning appearance is and then find all evening appearances
3756 * that take place before that */
3757 if (ipl == SE_MERCURY) {
3758 event_type = SE_HELIACAL_RISING;
3759 retflag = swe_heliacal_ut(dret[0], geopos, datm, dobs, obj_name, event_type, helflag, dret, serr);
3760 if (retflag == ERR) {
3761 do_printf(serr);
3762 return ERR;
3763 }
3764 tsave2 = dret[0];
3765 }
3766 repeat_mercury:
3767 /* evening first */
3768 event_type = SE_EVENING_FIRST;
3769 retflag = swe_heliacal_ut(tsave1, geopos, datm, dobs, obj_name, event_type, helflag, dret, serr);
3770 if (retflag == ERR) {
3771 do_printf(serr);
3772 return ERR;
3773 }
3774 if (ipl == SE_MERCURY && dret[0] > tsave2)
3775 continue;
3776 if (time_flag & (BIT_TIME_LMT | BIT_TIME_LAT)) {
3777 retc = ut_to_lmt_lat(dret[0], geopos, &(dret[0]), serr);
3778 if (retc != ERR) retc = ut_to_lmt_lat(dret[1], geopos, &(dret[1]), serr);
3779 if (retc != ERR) retc = ut_to_lmt_lat(dret[2], geopos, &(dret[2]), serr);
3780 if (retc == ERR) {
3781 do_printf(serr);
3782 return ERR;
3783 }
3784 }
3785 do_print_heliacal(dret, event_type, obj_name);
3786 }
3787 if (ipl == SE_MOON) {
3788 /* morning last */
3789 event_type = SE_MORNING_LAST;
3790 retflag = swe_heliacal_ut(dret[0], geopos, datm, dobs, obj_name, event_type, helflag, dret, serr);
3791 if (retflag == ERR) {
3792 do_printf(serr);
3793 return ERR;
3794 }
3795 if (time_flag & (BIT_TIME_LMT | BIT_TIME_LAT)) {
3796 retc = ut_to_lmt_lat(dret[0], geopos, &(dret[0]), serr);
3797 if (retc != ERR) retc = ut_to_lmt_lat(dret[1], geopos, &(dret[1]), serr);
3798 if (retc != ERR) retc = ut_to_lmt_lat(dret[2], geopos, &(dret[2]), serr);
3799 if (retc == ERR) {
3800 do_printf(serr);
3801 return ERR;
3802 }
3803 }
3804 do_print_heliacal(dret, event_type, obj_name);
3805 } else {
3806 /* heliacal setting (evening last) */
3807 event_type = SE_HELIACAL_SETTING;
3808 retflag = swe_heliacal_ut(dret[0], geopos, datm, dobs, obj_name, event_type, helflag, dret, serr);
3809 if (retflag == ERR) {
3810 do_printf(serr);
3811 return ERR;
3812 }
3813 if (time_flag & (BIT_TIME_LMT | BIT_TIME_LAT)) {
3814 retc = ut_to_lmt_lat(dret[0], geopos, &(dret[0]), serr);
3815 if (retc != ERR) retc = ut_to_lmt_lat(dret[1], geopos, &(dret[1]), serr);
3816 if (retc != ERR) retc = ut_to_lmt_lat(dret[2], geopos, &(dret[2]), serr);
3817 if (retc == ERR) {
3818 do_printf(serr);
3819 return ERR;
3820 }
3821 }
3822 do_print_heliacal(dret, event_type, obj_name);
3823 if ((0) && ipl == SE_MERCURY) {
3824 tsave1 = dret[0];
3825 goto repeat_mercury;
3826 }
3827 }
3828 }
3829 }
3830 return OK;
3831 }
3832
do_special_event(double tjd,int32 ipl,char * star,int32 special_event,int32 special_mode,double * geopos,double * datm,double * dobs,char * serr)3833 static int do_special_event(double tjd, int32 ipl, char *star, int32 special_event, int32 special_mode, double *geopos, double *datm, double *dobs, char *serr)
3834 {
3835 int retc = 0;
3836 /* risings, settings, meridian transits */
3837 if (special_event == SP_RISE_SET ||
3838 special_event == SP_MERIDIAN_TRANSIT)
3839 retc = call_rise_set(tjd, ipl, star, whicheph, geopos, serr);
3840 /* lunar eclipses */
3841 if (special_event == SP_LUNAR_ECLIPSE)
3842 retc = call_lunar_eclipse(tjd, whicheph, special_mode, geopos, serr);
3843 /* solar eclipses */
3844 if (special_event == SP_SOLAR_ECLIPSE)
3845 retc = call_solar_eclipse(tjd, whicheph, special_mode, geopos, serr);
3846 /* occultations by the moon */
3847 if (special_event == SP_OCCULTATION)
3848 retc = call_lunar_occultation(tjd, ipl, star, whicheph, special_mode, geopos, serr);
3849 /* heliacal event */
3850 if (special_event == SP_HELIACAL)
3851 retc = call_heliacal_event(tjd, ipl, star, whicheph, geopos, datm, dobs, serr);
3852 return retc;
3853 }
3854
hms_from_tjd(double tjd)3855 static char *hms_from_tjd(double tjd)
3856 {
3857 static char s[AS_MAXCH];
3858 double x;
3859 /* tjd may be negative, 0h corresponds to day number 9999999.5 */
3860 x = fmod(tjd, 1); /* may be negative ! */
3861 x = fmod(x + 1.5, 1); /* is positive day fraction */
3862 sprintf(s, "%s ", hms(x * 24, BIT_LZEROES));
3863 return s;
3864 }
3865
hms(double x,int32 iflag)3866 static char *hms(double x, int32 iflag)
3867 {
3868 static char s[AS_MAXCH], s2[AS_MAXCH], *sp;
3869 char *c = ODEGREE_STRING;
3870 x += 0.5 / 36000.0; /* round to 0.1 sec */
3871 strcpy(s, dms(x, iflag));
3872 sp = strstr(s, c);
3873 if (sp != NULL) {
3874 *sp = ':';
3875 strcpy(s2, sp + strlen(ODEGREE_STRING));
3876 strcpy(sp + 1, s2);
3877 *(sp + 3) = ':';
3878 *(sp + 8) = '\0';
3879 }
3880 return s;
3881 }
3882
do_printf(char * info)3883 static void do_printf(char *info)
3884 {
3885 #ifdef _WINDOWS
3886 fprintf(fp, info);
3887 #else
3888 fputs(info,stdout);
3889 #endif
3890 }
3891
3892 /* make_ephemeris_path().
3893 * ephemeris path includes
3894 * current working directory
3895 * + program directory
3896 * + default path from swephexp.h on current drive
3897 * + on program drive
3898 * + on drive C:
3899 */
make_ephemeris_path(char * argv0,char * path)3900 static int make_ephemeris_path(char *argv0, char *path)
3901 {
3902 char *sp;
3903 char *dirglue = DIR_GLUE;
3904 size_t pathlen = 0;
3905 /* current working directory */
3906 sprintf(path, ".%c", *PATH_SEPARATOR);
3907 /* program directory */
3908 sp = strrchr(argv0, *dirglue);
3909 if (sp != NULL) {
3910 pathlen = sp - argv0;
3911 if (strlen(path) + pathlen < AS_MAXCH-2) {
3912 strncat(path, argv0, pathlen);
3913 sprintf(path + strlen(path), "%c", *PATH_SEPARATOR);
3914 }
3915 }
3916 #if MSDOS
3917 {
3918 char *cpos[20];
3919 char s[2 * AS_MAXCH], *s1 = s + AS_MAXCH;
3920 char *sp[3];
3921 int i, j, np;
3922 strcpy(s1, SE_EPHE_PATH);
3923 np = cut_str_any(s1, PATH_SEPARATOR, cpos, 20);
3924 /*
3925 * default path from swephexp.h
3926 * - current drive
3927 * - program drive
3928 * - drive C
3929 */
3930 *s = '\0';
3931 /* current working drive */
3932 sp[0] = getcwd(NULL, 0);
3933 if (sp[0] == NULL) {
3934 printf("error in getcwd()\n");
3935 exit(1);
3936 }
3937 if (*sp[0] == 'C')
3938 sp[0] = NULL;
3939 /* program drive */
3940 if (*argv0 != 'C' && (sp[0] == NULL || *sp[0] != *argv0))
3941 sp[1] = argv0;
3942 else
3943 sp[1] = NULL;
3944 /* drive C */
3945 sp[2] = "C";
3946 for (i = 0; i < np; i++) {
3947 strcpy(s, cpos[i]);
3948 if (*s == '.') /* current directory */
3949 continue;
3950 if (s[1] == ':') /* drive already there */
3951 continue;
3952 for (j = 0; j < 3; j++) {
3953 if (sp[j] != NULL && strlen(path) + 2 + strlen(s) < AS_MAXCH-1) {
3954 sprintf(path + strlen(path), "%c:%s%c", *sp[j], s, *PATH_SEPARATOR);
3955 }
3956 }
3957 }
3958 }
3959 #else
3960 if (strlen(path) + strlen(SE_EPHE_PATH) < AS_MAXCH-1)
3961 strcat(path, SE_EPHE_PATH);
3962 #endif
3963 return OK;
3964 }
3965
remove_whitespace(char * s)3966 static void remove_whitespace(char *s)
3967 {
3968 char *sp, s1[AS_MAXCH];
3969 if (s == NULL) return;
3970 for (sp = s; *sp == ' '; sp++)
3971 ;
3972 strcpy(s1, sp);
3973 while (*(sp = s1 + strlen(s1) - 1) == ' ')
3974 *sp = '\0';
3975 strcpy(s, s1);
3976 }
3977
jd_to_time_string(double jut,char * stimeout)3978 static void jd_to_time_string(double jut, char *stimeout)
3979 {
3980 double t2;
3981 t2 = jut + 0.5 / 3600000.0; // rounding to millisec
3982 sprintf(stimeout, " % 2d:", (int) t2); // hour
3983 t2 = (t2 - (int32) t2) * 60;
3984 sprintf(stimeout + strlen(stimeout), "%02d:", (int) t2); // min
3985 t2 = (t2 - (int32) t2) * 60;
3986 sprintf(stimeout + strlen(stimeout), "%02d", (int) t2); // sec
3987 t2 = (t2 - (int32) t2) * 1000;
3988 if ((int32) t2 > 0) {
3989 sprintf(stimeout + strlen(stimeout), ".%03d", (int) t2); // millisec, if > 0
3990 }
3991 }
3992
3993 #if MSDOS
3994 /**************************************************************
3995 cut the string s at any char in cutlist; put pointers to partial strings
3996 into cpos[0..n-1], return number of partial strings;
3997 if less than nmax fields are found, the first empty pointer is
3998 set to NULL.
3999 More than one character of cutlist in direct sequence count as one
4000 separator only! cut_str_any("word,,,word2",","..) cuts only two parts,
4001 cpos[0] = "word" and cpos[1] = "word2".
4002 If more than nmax fields are found, nmax is returned and the
4003 last field nmax-1 rmains un-cut.
4004 **************************************************************/
cut_str_any(char * s,char * cutlist,char * cpos[],int nmax)4005 static int cut_str_any(char *s, char *cutlist, char *cpos[], int nmax)
4006 {
4007 int n = 1;
4008 cpos [0] = s;
4009 while (*s != '\0') {
4010 if ((strchr(cutlist, (int) *s) != NULL) && n < nmax) {
4011 *s = '\0';
4012 while (*(s + 1) != '\0' && strchr (cutlist, (int) *(s + 1)) != NULL) s++;
4013 cpos[n++] = s + 1;
4014 }
4015 if (*s == '\n' || *s == '\r') { /* treat nl or cr like end of string */
4016 *s = '\0';
4017 break;
4018 }
4019 s++;
4020 }
4021 if (n < nmax) cpos[n] = NULL;
4022 return (n);
4023 } /* cutstr */
4024 #endif
4025
our_strcpy(char * to,char * from)4026 static char *our_strcpy(char *to, char *from)
4027 {
4028 char *sp, s[AS_MAXCH];
4029 if (*from == '\0') {
4030 *to = '\0';
4031 return to;
4032 }
4033 if (strlen(from) < AS_MAXCH) {
4034 strcpy(s, from);
4035 strcpy(to, s);
4036 } else {
4037 sp = strdup(from);
4038 if (sp == NULL) {
4039 strcpy(to, from);
4040 } else {
4041 strcpy(to, sp);
4042 free(sp);
4043 }
4044 }
4045 return to;
4046 }
4047