source: Sophya/trunk/SophyaLib/NTools/datime.c@ 629

Last change on this file since 629 was 517, checked in by ansari, 26 years ago

no message

File size: 34.0 KB
Line 
1#include "machdefs.h"
2#include <stdlib.h>
3#include <stdio.h>
4#include <string.h>
5#include <ctype.h>
6#include <math.h>
7#ifdef __MWERKS__
8 #include "mwerksmath.h"
9// #include "unixmac.h"
10#endif
11
12#define DATIMEPRIVEE
13#include "strutil.h"
14#include "nbconst.h"
15#include "datime.h"
16
17
18/* Fonctions de calcule de date et temps (Heure) */
19/*
20++
21 Module Dates (C)
22 Lib LibsUtil
23 include datime.h
24
25 Ce groupe de fonctions permettent de manipuler des dates et heures.
26 En particulier, il est possible de calculer l'ecart (en nombre de jours
27 ou de secondes separant deux dates, ou le temps sideral correspondant
28 a une date et heure legale. Deux structures simples sont definies
29 afin de faciliter le passage des arguments entre differentes fonctions:
30 - *JMA* : Jour, Mois, Annee
31 - *HMS* : Heure, Minutes, Secondes
32
33--
34*/
35/*
36++
37 Links Voir aussi:
38 Temps Sideral, Universel (C)
39--
40*/
41/*
42++
43Titre Quelques Macros
44--
45*/
46
47/*
48++
49 StrgtoJMA(strg, jma)
50 Decodage d'une chaine de caracteres "strg" sous forme de "12/3/97" en structure
51 JMA "jma".
52 JMAtoStrg(jma, strg)
53 Ecrit le contenu de la structure JMA "jma" sous forme de "jj/mm/aaaa" dans la chaine
54 de caracteres "strg"
55 JMAtoStrgLong(jma, strg)
56 Ecriture de date en format long , par ex : "Jeudi , 19 Juin 1997"
57 StrgtoHMS(strg, hms)
58 Decodage d'une chaine de caracteres "strg" sous forme de "hh:mm:ss 10:43:60.5"
59 en structure HMS "hms".
60 HMStoStrg(hms, strg)
61 Ecrit le contenu de la structure HMS "hms" sous forme de "hh:mm:ss" dans la chaine
62 de caracteres "strg"
63 StrtoHMS(strg, hms)
64 Decodage d'une chaine de caracteres "strg" sous forme de "hh:mm:ss 10:43:60.5"
65 en structure HMS "hms" avec gestion des signes (decodage correcte de -44:30:05.23)
66 HMStoStr(hms, strg)
67 Ecrit le contenu de la structure HMS "hms" sous forme de "hh:mm:ss" dans la chaine
68 de caracteres "strg" avec gestion des signes
69--
70*/
71
72void StrtoHMS(char *s,HMS* h)
73/* On ne peut pas ecrire 1:-05:-45 pas gere et debile! (mn et sec >=0.) cmv 12/8/97 */
74{
75int imoins, i2pt;
76h->Heures = h->Minutes = 0; h->Secondes = 0.;
77/* decodage de la chaine de characteres */
78sscanf(s,"%d:%d:%lf", &(h->Heures),&(h->Minutes),&(h->Secondes));
79/* cas du type 5:34:45.5 */
80if( h->Heures > 0 ) return;
81/* cas du type -5:34:45.5 */
82if( h->Heures < 0 ) {
83 h->Minutes *= -1; h->Secondes *= -1.;
84 return;
85}
86/* Ici les Heures sont donc nulles "0:..." ou "-0:..." */
87/* on n'a donne que l'heure nulle cas "0" ou "-0" sans ":" (idiot mais bon!) */
88i2pt = posc(s,':');
89if( i2pt < 0 ) return;
90/* Ici on a au moins un ou 2 ":" */
91/* cas du type 0:34:45.5 ou 0:0:45.5 sans signe - */
92imoins = posc(s,'-');
93if( imoins < 0 ) return;
94/* cas du type -0:34:45.5 ou -0:0:45.5, le signe - n'est pas decode par sscanf */
95if( imoins < i2pt ) { /* on a donne -0:34 ou -0:34:45.5 */
96 h->Minutes *= -1; h->Secondes *= -1.;
97 return;
98}
99/* Ici il reste les cas mal codes ou le signe "-" est n'importe ou! */
100/* On laisse le decodage par sscanf */
101return;
102}
103
104
105void HMStoStr(HMS h,char *s)
106/* L'ecriture est forcee a h:mn:sec avec mn et sec >=0 cmv 12/8/97 */
107{
108double v;
109v = HMStoH(h);
110h = DoubletoHMS(v);
111if( v < 0. ) {
112 h.Heures *= -1; h.Minutes *= -1; h.Secondes *= -1.;
113 sprintf(s,"-%02d:%02d:%04.1f", h.Heures,h.Minutes,h.Secondes);
114} else {
115 sprintf(s,"%02d:%02d:%04.1f", h.Heures,h.Minutes,h.Secondes);
116}
117}
118
119/* Nouvelle-Fonction */
120int NbJourMois(int a, int m)
121/* Retourne le nombre de jours dans le mois m, annee a */
122{
123if(a<100) a+=1900;
124if( m<1 || m>12 ) return(-1);
125m--;
126if(m!=1) return(NbJoMois[m]);
127/* traitement du mois de fevrier bissextile ou non */
128if( ((a%4==0) && (a%100!=0))
129 || (a%400==0) ) return(NbJoMois[m]+1);
130else return(NbJoMois[m]);
131}
132
133/* Nouvelle-Fonction */
134long JMAtoJ(JMA jma)
135/* Calcule le Nb. de jours ecoules depuis 0 Jan 1901 */
136/* Si annee < 100 On considere annee = annee+1900 (1900-1999) */
137
138{
139long rc,nban;
140int i;
141
142/* Protection */
143if ( (jma.Mois < 1) || (jma.Mois > 12) ) jma.Mois = 1;
144if ( (jma.Jour < 1) || (jma.Jour > 31) ) jma.Jour = 1;
145
146if (jma.Annee < 100) jma.Annee+=1900;
147nban = jma.Annee - 1901;
148
149/* Annee Normale = 365 Jours - Bissextile = 366 Jours */
150/* 366 Jours : Divisible par 4 et non divisible par 100 (1972,1976) */
151/* : Divisible par 400 (1900 Non-Biss. et 2000 Biss.) */
152
153if (nban >= 0)
154 rc = nban*365 + (nban/4) - (nban/100) + ((nban+300)/400);
155else
156 rc = nban*365 + (nban/4) - (nban/100) + ((nban-300)/400);
157
158/* Comptabilisation du nb. des jours des mois ecoule */
159for (i=1; i<jma.Mois; i++) rc += NbJourMois(jma.Annee, i);
160rc += jma.Jour;
161
162return(rc);
163}
164
165/* Nouvelle-Fonction */
166JMA JtoJMA(long j)
167/* retourne la date correspondant a j jours ecoules depuis le 0/1/1901 */
168{
169long i;
170int inc,m;
171JMA buf;
172
173inc = (j>=0) ? 1 : -1;
174buf.Jour = buf.Mois = 1; buf.Annee = 1901;
175
176/* recherche de l'annee */
177for(;;) {
178 i = j-JMAtoJ(buf);
179 if( i==0 ) return(buf);
180 if( inc*i<0 ) {
181 /* printf("...annee j=%ld scan=%ld test=%ld\n"
182 ,j,JMAtoJ(buf),inc*i); */
183 break;
184 }
185 buf.Annee += inc;
186}
187if(inc>0) buf.Annee--;
188/* printf("Annee trouvee: %d\n",buf.Annee); */
189
190/* recherche du mois */
191for(m=1;m<=12;m++) {
192 buf.Mois=m;
193 i = j-JMAtoJ(buf);
194 if( i==0 ) return(buf);
195 if(i<0) {
196 /* printf("...mois j=%ld scan=%ld test=%ld\n"
197 ,j,JMAtoJ(buf),i); */
198 buf.Mois--;
199 break;
200 }
201}
202/* printf("Mois trouve: %d\n",buf.Mois); */
203
204/* recherche du jour */
205for(m=1;m<=31;m++) {
206 buf.Jour=m;
207 i = j-JMAtoJ(buf);
208 if( i==0 ) return(buf);
209}
210/* printf("Jour trouve: %d\n",buf.Jour); */
211
212/* on a pas trouve !!! normallement impossible! */
213printf("**** JtoJMA_Erreur: date non trouvee j=%ld ****\n",j);
214buf.Jour = buf.Mois = 1; buf.Annee = 0;
215return(buf);
216}
217
218/* Nouvelle-Fonction */
219int NumJour(JMA jma)
220/* Cette fonction calcule le numero de jour de la semaine */
221/* (Lundi=1 .. Dimanche=7) */
222
223{
224int l;
225
226/* On appelle JMAtoJ() pour avoir le nb de jours ecoule depuis 1/01/1901 */
227/* Le premier janvier 1901 etait un Mardi (JMAtoJ = 1 , NumJour = 2) */
228l = JMAtoJ(jma)%7 ;
229if (l >= 0) return(l+1);
230else return(l+8);
231}
232
233/* Nouvelle-Fonction */
234JMA JsmmaatoJMA(int nj, int mm, int aa)
235
236/* Pour "nj=kj*10+jj", calcule la date correspondant au "kj" eme
237 jour "jj" (1 .. 7 - Lundi .. Dimanche) du mois mm , annee aa
238 Si "nj" est negatif -> Dernier jour "jj". Si "kj=0", premier jour "jj".
239 Si "kj" ne peut etre atteint -> Dernier jour "jj".
240 Si annee "aa" <= 99 -> annee de 1901 - 1999 */
241
242{
243JMA jma;
244int i,k;
245int jj, kj;
246
247kj = nj/10;
248jj = nj%10;
249if (kj < 0) kj = 0;
250
251if (mm > 12) mm = 12;
252if (mm < 1 ) mm = 1;
253if (aa < 100) aa += 1900;
254jma.Mois = mm;
255jma.Annee = aa;
256
257if (kj == 0) {
258 if (jj >= 0) { /* Premier jour du mois */
259 if (jj > 7) jj = 7;
260 if (jj < 1) jj = 1;
261 jma.Jour = 1;
262 i = NumJour(jma);
263 if (i > jj) jma.Jour = (7-i)+jj+1;
264 else jma.Jour = 1+(jj-i);
265 }
266 else { /* Dernier jour du mois */
267 jj = -jj;
268 if (jj > 7) jj = 7;
269 if (jj < 1) jj = 1;
270 jma.Jour = NbJourMois(aa,mm);
271 i = NumJour(jma);
272 if (i >= jj) jma.Jour -= (i-jj);
273 else jma.Jour -= (7+i-jj);
274 }
275}
276else { /* On demande le kj eme jour jj du mois */
277 if (jj > 7) jj = 7;
278 if (jj < 1) jj = 1;
279 k = 0;
280 for(i=1; i<=NbJourMois(aa,mm); i++) {
281 jma.Jour = i;
282 if (NumJour(jma) == jj) k++;
283 if (k == kj) break;
284 }
285/* On n'a pas kj jour jj ds ce mois, on renvoie le dernier jj */
286 if (k<kj) return(JsmmaatoJMA(-jj, mm, aa));
287}
288
289return(jma);
290}
291
292JMA JApmmaatoJMA(int nj, int mm, int aa)
293/* Pour "nj=dd*10+jj", calcule la date correspondant au jour
294 "jj" (1 .. 7 - Lundi .. Dimanche) du mois "mm" , annee "aa"
295 se trouvant apres (ou egal a) la date "dd/mm/aa" .
296 Si annee "aa" <= 99 -> annee de 1901 - 1999 */
297{
298int i, dd, jj;
299JMA jma;
300
301if (nj < 0) nj = -nj;
302if (nj < 1) nj = 1;
303dd = nj/10;
304jj = nj%10;
305for(i=1; i<=6; i++) {
306 jma = JsmmaatoJMA(jj+10*i, mm, aa);
307 if(jma.Jour >= dd) return(jma);
308}
309printf("JApmmaatoJMA(%d...) Erreur / Pas de jour %d >= %d/%d/%d\n", nj, jj,dd, mm, aa);
310jma.Jour = jma.Mois = jma.Annee = 0;
311return(jma);
312}
313
314/* Nouvelle-Fonction */
315double HMStoH(HMS hms)
316
317/* Cette fonction calcule le nb d'heures en decimales pour hms */
318/* Heures/Minutes/Secondes peuvent etre +/- avec toutes les combi possibles */
319{
320return ((double)hms.Heures + ((double)hms.Minutes)/60. + hms.Secondes/3600.);
321}
322
323
324
325/* Nouvelle-Fonction */
326double HMStoSec(HMS hms)
327/* Calcul du nb de secondes pour hms */
328/* Heures/Minutes/Secondes peuvent etre +/- avec toutes les combi possibles */
329{
330return ((double)hms.Heures*3600. + (double)hms.Minutes*60. + hms.Secondes);
331}
332
333/* Nouvelle-Fonction */
334long DatetoSec(char const* date, char const* heure)
335/* Calcul du nombre de secondes correspondant a date/heure */
336/* a partir de l'origine 1er Janv 1990 0H00 */
337/* Date entre 1930 - 2050 */
338{
339JMA jma,jma0;
340HMS hms;
341long dj;
342
343StrgtoJMA(date,jma); StrgtoHMS(heure,hms);
344StrgtoJMA("01/01/1990",jma0);
345dj = JMAtoJ(jma) - JMAtoJ(jma0);
346if (dj < -24000) return(-0x7FFFFFFF);
347if (dj > 24000) return(0x7FFFFFFF);
348return( dj*3600*24 + (int)(HMStoSec(hms)) );
349}
350
351/* Nouvelle-Fonction */
352long DatetoSecOff(char const* date, char const* heure)
353{
354JMA jma,jma0;
355HMS hms, hms0;
356long dj,tdec;
357
358StrgtoJMA(date,jma); StrgtoHMS(heure,hms);
359StrgtoJMA("01/01/1990",jma0); StrgtoHMS("00:00:00",hms0);
360dj = JMAtoJ(jma) - JMAtoJ(jma0);
361if (dj < -24000) return(-0x7FFFFFFF);
362if (dj > 24000) return(0x7FFFFFFF);
363
364/* H_vrai+TOff=H_leg : H_vrai-H_vrai0 = (H_leg-TOff)-(H_leg0-TOff0) */
365tdec = 3600*(TLegOffset(jma, hms)-TLegOffset(jma0, hms0));
366
367dj = dj*3600*24 + (int)(HMStoSec(hms)) - tdec;
368return(dj);
369}
370
371/* Nouvelle-Fonction */
372HMS HtoHMS(double h)
373/* Heures decimales en H , M , Sec Il y a un modulo 24 */
374/* la sortie est telle que 0 <= heures < 24 */
375{
376HMS hms;
377
378while (h < 0.) h += 24.;
379
380hms.Heures = (int) h;
381h = (h-hms.Heures)*60.;
382hms.Minutes = (int) h;
383hms.Secondes = (h-hms.Minutes)*60.;
384hms.Heures = hms.Heures % 24;
385
386return(hms);
387}
388
389/* Nouvelle-Fonction */
390HMS DtoDMS(double h)
391/* Degres decimaux en D , M , Sec Il y a un modulo 360 */
392/* la sortie est telle que -180 < deg <= 180. */
393{
394HMS hms;
395int sgn = 1;
396
397while (h > 180.) h -= 360.;
398while (h <= -180.) h += 360.;
399if( h < 0. ) {h *= -1.; sgn = -1;}
400
401hms.Heures = (int) h;
402h = (h-hms.Heures)*60.;
403hms.Minutes = (int) h;
404hms.Secondes = (h-hms.Minutes)*60.;
405hms.Heures = hms.Heures % 360;
406
407if( sgn<0 ) {hms.Heures *= -1; hms.Minutes *= -1; hms.Secondes *= -1.;}
408return(hms);
409}
410
411/* Nouvelle-Fonction */
412HMS DoubletoHMS(double h)
413/* Degres/Heures decimaux en D/H , M , Sec SANS modulo 360/24 */
414{
415HMS hms;
416
417int sgn = 1;
418if( h < 0. ) {h *= -1.; sgn = -1;}
419
420hms.Heures = (int) h;
421h = (h-hms.Heures)*60.;
422hms.Minutes = (int) h;
423hms.Secondes = (h-hms.Minutes)*60.;
424hms.Heures = hms.Heures;
425
426if( sgn<0 ) {hms.Heures *= -1; hms.Minutes *= -1; hms.Secondes *= -1.;}
427return(hms);
428}
429
430/*
431++
432 Titre Calcule sur dates et heures
433--
434*/
435
436/*
437++
438int NbJourMois(int a, int m)
439 Retourne le nombre de jours dans le mois "m", annee "a".
440long JMAtoJ(JMA jma)
441 Calcule le Nb. de jours ecoules depuis 0 Jan 1901,
442 Si annee < 100 On considere annee = annee+1900 (1900-1999)
443JMA JtoJMA(long j)
444 Retourne la date correspondant a un nombre de jours "j"
445 ecoules depuis le 0 Jan 1901.
446int NumJour(JMA jma)
447 Calcule le numero du jour de la semaine correspondant a la date indiquee
448 par la structure JMA "jma". (Lundi=1 .. Dimanche=7)
449JMA JsmmaatoJMA(int nj, int mm, int aa)
450 Pour "nj=kj*10+jj", calcule la date correspondant au "kj" eme
451 jour "jj" (1 .. 7 - Lundi .. Dimanche) du mois mm , annee aa
452 Si "nj" est negatif -> Dernier jour "jj". Si "kj=0", premier jour "jj".
453 Si "kj" ne peut etre atteint -> Dernier jour "jj".
454 Si annee "aa" <= 99 -> annee de 1901 - 1999
455JMA JApmmaatoJMA(int nj, int mm, int aa)
456 Pour "nj=dd*10+jj", calcule la date correspondant au jour
457 "jj" (1 .. 7 - Lundi .. Dimanche) du mois "mm" , annee "aa"
458 se trouvant apres (ou egal a) la date "dd/mm/aa" .
459 Si annee "aa" <= 99 -> annee de 1901 - 1999
460double HMStoH(HMS hms)
461 Conversion en heures decimales de la structure "hms"
462double HMStoSec(HMS hms)
463 Conversion en secondes de la structure HMS "hms"
464HMS HtoHMS(double h)
465 Conversion heures en decimales en structure HMS
466HMS DtoDMS(double h)
467 Conversion degres decimaux en structure HMS (deg,min,sec)
468HMS DoubletoHMS(double h)
469 Conversion degres/heures decimaux en structure HMS (deg/heure,min,sec)
470 SANS modulo 360/24
471long DatetoSec(char const* date, char const* heure)
472 Calcule le nombre de secondes correspondant a "date,heure"
473 a partir de l'origine 1er Janv 1990 0H00 ( Date entre 1930 - 2050 )
474long DatetoSecOff(char const* date, char const* heure)
475 Calcule le nombre de secondes correspondant a "date,heure"
476 a partir de l'origine 1er Janv 1990 0H00 ( Date entre 1930 - 2050 )
477 en tenant compte des decalages des heures legales.
478--
479*/
480
481
482
483/* Variables contenant les informations pour conversion */
484/* Temps Legal <> TU <> Temps solaire moyen */
485
486static double TSolMOff = 0.; /* TSol.Moyen - TU (en heures) */
487static char SiteName[128]; /* Nom de site */
488
489static struct
490 {
491 char fg; /* Flag (=A ou B) Type de specification */
492 int jj[2]; /* Numero ou date du jour */
493 int mm[2]; /* Mois */
494 int Off[3]; /* Decalage = TLegal - TU */
495 } TLegInfo = {'\0',{0,0},{0,0},{0,0,0}} ;
496
497static HMS HChgH = {0,0,0.} ; /* Heure de changement d'heure */
498
499/*
500++
501 Module Temps Sideral, Universel (C)
502 Lib LibsUtil
503 include datime.h
504
505 Ce groupe de fonctions permet convertir differentes heures:
506 les temps local, le temps universel (TU), le temps legal (TLeg). Il est
507 aussi possible de calculer le temps sideral (TS).
508--
509*/
510/*
511++
512 Links Voir aussi:
513 Dates (C)
514--
515*/
516/*
517++
518 Titre T-Legal, TU
519 Conversion entre l'heure legale, le temps universel (TU) et le temps solaire
520 moyen
521--
522*/
523
524/*
525++
526SetTSolMOff(double tsmoff, char *name)
527 Definit pour un site donne, le decalage T Sol. Moyen - TU (en heures) "tmsoff"
528 et le nom du site
529PrtTSolMOff()
530 Imprime les infos definis par "SetTSolMOff()"
531--
532*/
533
534/* Nouvelle-Fonction */
535void SetTSolMOff(double tsmoff, char *name)
536
537/* Offset T Sol. Moyen - TU (en hh:mm:ss) , Nom de site name */
538{
539TSolMOff = tsmoff;
540strcpy(SiteName,name);
541return;
542}
543
544
545/* Nouvelle-Fonction */
546void PrtTSolMOff(void)
547{
548printf(" **Site : %s , T Sol.Moy.-TU %g (heures) \n",SiteName,TSolMOff);
549return;
550}
551
552/*
553++
554int SetTLegOff(char* dtz, char* dhz)
555 Initialise la structure TLegInfo. L'argument *dtz peut etre
556 sous 2 forme differentes et permet de definir 3 plages
557 de temps avec 3 decalages horaires (nb entier).
558 Les dates pivots sont les jours de changement d'heure a
559 l'heure indiquee par dhz ("hh:mm:ss"), si "dhz==NULL ou
560 chaine vide, alors changement a 0 heure.
561| 1) dtz = "A Off1 j1/m1 Off2 j2/m2 Off3"
562| jx/mx Jour/Mois - Offx : Decalage TLeg - TU
563| pour 01/01 <= jj/mm < j1/m1 -> Off1
564| j1/m1 <= jj/mm < j2/m2 -> Off2
565| j2/m2 <= jj/mm <= 31/12 -> off3
566| 2) dtz = "B off1 Js1/m1 Off2 Js2/m2 Off3"
567| Jsx/mx Numero du jour de la semaine (Js) et numero du mois
568| Js = 1 .. 7 (Lundi .. Dimanche) : 1er Jour Js du mois
569| = -1 .. -7 Dernier jour Js du mois
570| = k1 .. k7 k eme jour Js du mois
571| 3) dtz = "C off1 Js1/m1 Off2 Js2/m2 Off3"
572| Jsx/mx Numero du jour de la semaine (Js) et numero du mois avec une
573| precision d'une date butee inferieurs
574| Js = ddj .. (j=1..7 Lundi .. Dimanche) : 1er Jour j du mois
575| apres (ou egal a) la date dd/mm
576 Retour : 0 si Ok -1 Sinon
577void PrtTLegInfo(int aa)
578 Imprime les informations de decalage TLeg , TU - Annee "aa"
579int TLegOffset(JMA date, HMS heure)
580 Indique le decalage (en nombre d'heures) TLeg-TU
581--
582*/
583
584/* Nouvelle-Fonction */
585int SetTLegOff(char* dtz, char* dhz)
586
587
588{
589int is;
590
591HChgH.Heures = HChgH.Minutes = 0;
592HChgH.Secondes = 0.;
593
594if (dhz) StrgtoHMS(dhz, HChgH);
595if ( (HChgH.Heures < 0) || (HChgH.Heures > 23) )
596 { HChgH.Heures = HChgH.Minutes = 0; HChgH.Secondes = 0.; }
597
598is = sscanf(dtz,"%c %d %d/%d %d %d/%d %d", &(TLegInfo.fg), &(TLegInfo.Off[0]), &(TLegInfo.jj[0]), &(TLegInfo.mm[0]), &(TLegInfo.Off[1]), &(TLegInfo.jj[1]), &(TLegInfo.mm[1]), &(TLegInfo.Off[2]) );
599
600if ( (TLegInfo.fg != 'A') && (TLegInfo.fg != 'B')
601 && (TLegInfo.fg != 'C') ) TLegInfo.fg = '\0';
602if (is != 8) TLegInfo.fg = '\0';
603if (TLegInfo.fg != '\0') return(0);
604else { printf("SetTLegOff()/Erreur %s \n", dtz); return(-1); }
605}
606
607/* Nouvelle-Fonction */
608void PrtTLegInfo(int aa)
609
610/* Cette fonction print les infos de TLegInfo */
611{
612int i,k,j;
613char ForL[2][128];
614JMA jma;
615
616if (TLegInfo.fg == '\0') printf (" !!! TlegInfo not set or Error \n");
617
618
619HMStoStrg(HChgH, ForL[0]);
620printf(" ==== TLegInfo : Changement d'heure a %s \n", ForL[0]);
621ForL[0][0] = '\0';
622if (TLegInfo.fg == 'A')
623 {
624 printf(" ---- TLegInfo : Selection par date (Flag = A) ---- \n");
625 printf(" -> 01 Janvier - %d %s : %d heures \n",TLegInfo.jj[0], NomMo[TLegInfo.mm[0]-1], TLegInfo.Off[0]);
626 printf(" -> %d %s - %d %s : %d heures \n",TLegInfo.jj[0], NomMo[TLegInfo.mm[0]-1],TLegInfo.jj[1], NomMo[TLegInfo.mm[1]-1],TLegInfo.Off[1]);
627 printf(" -> %d %s - 31 Decembre : %d heures \n",TLegInfo.jj[1], NomMo[TLegInfo.mm[1]-1], TLegInfo.Off[2]);
628 }
629
630if (TLegInfo.fg == 'B')
631 {
632 printf(" --- TLegInfo : Selection par Jour de semaine et mois (Flag = B) --- \n");
633 for (i=0;i<2;i++) {
634 k = TLegInfo.jj[i]/10;
635 if (k < 1) k = 1;
636 j = TLegInfo.jj[i]%10;
637 if (j < 0) j = -j;
638 if (j < 1) j = 1;
639 if (j > 7) j = 7;
640 if (TLegInfo.jj[i] < 0) sprintf(ForL[i]," Dernier %s de %s (KJ=%d)", NomJo[j-1],
641 NomMo[TLegInfo.mm[i]-1], TLegInfo.jj[i]);
642 else sprintf(ForL[i]," %d eme %s de %s (KJ=%d)", k, NomJo[j-1],
643 NomMo[TLegInfo.mm[i]-1], TLegInfo.jj[i]);
644 }
645 printf(" Butees: %s - %s \n", ForL[0], ForL[1]);
646 for (i=0;i<2;i++) {
647 jma = JsmmaatoJMA(TLegInfo.jj[i], TLegInfo.mm[i], aa);
648 JMAtoStrg(jma, ForL[i]);
649 }
650 printf(" -> 01 Janvier %d - %s : %d heures \n",aa, ForL[0], TLegInfo.Off[0]);
651 printf(" -> %s - %s : %d heures \n",ForL[0], ForL[1], TLegInfo.Off[1]);
652 printf(" -> %s - 31 Decembre %d : %d heures \n",ForL[1], aa, TLegInfo.Off[2]);
653 }
654
655if (TLegInfo.fg == 'C')
656 {
657 printf(" --- TLegInfo : Selection par Jour de semaine >= dd et mois (Flag = C) --- \n");
658 for (i=0;i<2;i++) {
659 k = TLegInfo.jj[i]/10;
660 j = TLegInfo.jj[i]%10;
661 if (j < 0) j = -j;
662 if (j < 1) j = 1;
663 if (j > 7) j = 7;
664 sprintf(ForL[i],"1er %s apres le %d %s ", NomJo[j-1], k, NomMo[TLegInfo.mm[i]-1]);
665 }
666 printf(" Butees: %s - %s \n", ForL[0], ForL[1]);
667 for (i=0;i<2;i++) {
668 jma = JApmmaatoJMA(TLegInfo.jj[i], TLegInfo.mm[i], aa);
669 JMAtoStrg(jma, ForL[i]);
670 }
671 printf(" -> 01 Janvier %d - %s : %d heures \n",aa, ForL[0], TLegInfo.Off[0]);
672 printf(" -> %s - %s : %d heures \n",ForL[0], ForL[1], TLegInfo.Off[1]);
673 printf(" -> %s - 31 Decembre %d : %d heures \n",ForL[1], aa, TLegInfo.Off[2]);
674 }
675
676return;
677}
678
679
680/* Nouvelle-Fonction */
681int TLegOffset(JMA date, HMS heure)
682
683/* Cette fonction calcule l'offset TLeg - TU pour la date donnee date */
684{
685JMA jma1,jma2;
686long nbj1,nbj2,nbj;
687double decj1, decj2, decj;
688int Off;
689
690/* On recherche l'offset correspondant a la date */
691jma1.Annee = jma2.Annee = date.Annee;
692if (TLegInfo.fg == 'A')
693/* Periode defini par des dates ds les mois */
694 {
695 jma1.Mois = TLegInfo.mm[0];
696 jma1.Jour = TLegInfo.jj[0];
697 jma2.Mois = TLegInfo.mm[1];
698 jma2.Jour = TLegInfo.jj[1];
699 }
700if (TLegInfo.fg == 'B')
701/* Periodes definies par jour de semaine et mois */
702 {
703 jma1 = JsmmaatoJMA(TLegInfo.jj[0],TLegInfo.mm[0],date.Annee);
704 jma2 = JsmmaatoJMA(TLegInfo.jj[1],TLegInfo.mm[1],date.Annee);
705 }
706
707if (TLegInfo.fg == 'C')
708/* Periodes definies par jour de semaine et mois avec date inferieure */
709 {
710 jma1 = JApmmaatoJMA(TLegInfo.jj[0],TLegInfo.mm[0],date.Annee);
711 jma2 = JApmmaatoJMA(TLegInfo.jj[1],TLegInfo.mm[1],date.Annee);
712 }
713
714nbj = JMAtoJ(date);
715nbj1 = JMAtoJ(jma1);
716nbj2 = JMAtoJ(jma2);
717
718/* Les dates pivots sont les jours de changement d'heure a Heure= HChgH */
719decj = (double)nbj + HMStoH(heure)/24.;
720decj1 = (double)nbj1 + HMStoH(HChgH)/24.;
721decj2 = (double)nbj2 + HMStoH(HChgH)/24.;
722
723if ( decj < decj1 ) Off = TLegInfo.Off[0];
724else if (decj < decj2) Off = TLegInfo.Off[1];
725 else Off = TLegInfo.Off[2];
726
727/* Impression de debugging */
728/*Debug-- printf (" @@@ Debug_TLegOffset (date= %d/%d/%d) : %d H. \n",date.Jour,date.Mois,date.Annee,Off);
729printf(" @-1 %s , %d %s %d \n", NomJo[NumJour(jma1)-1], jma1.Jour, NomMo[jma1.Mois-1], jma1.Annee) ;
730printf(" @-2 %s , %d %s %d \n", NomJo[NumJour(jma2)-1], jma2.Jour, NomMo[jma2.Mois-1], jma2.Annee) ; --*/
731
732
733return (Off);
734}
735
736/*
737++
738HMS TLegtoTU(JMA date, HMS TLeg)
739 Conversion TLegal en TU pour une date donnée
740HMS TUtoTLeg(JMA date, HMS TLeg)
741 Conversion TU en TLegal pour une date donnée. L'heure légale calculée
742 peut être erronnée autour des dates de changement d'heure.
743HMS TUtoTSolM(HMS TU)
744 Conversion TU en temps solaire moyen
745--
746*/
747
748/* Nouvelle-Fonction */
749HMS TLegtoTU(JMA date, HMS hmsz)
750
751/* Cette fonction transforme le temps legal en temps universel */
752{
753HMS TU;
754int Off;
755
756Off = TLegOffset(date, hmsz);
757TU = hmsz;
758TU.Heures -= Off;
759if (TU.Heures >= 24) TU.Heures -= 24;
760if (TU.Heures < 0) TU.Heures += 24;
761return(TU);
762}
763
764
765
766/* Nouvelle-Fonction */
767HMS TUtoTLeg(JMA date, HMS hms)
768
769/* Cette fonction transforme le temps universel en temps legal */
770{
771HMS TLeg;
772int Off;
773
774Off = TLegOffset(date, hms); /* $CHECK$ */
775TLeg = hms;
776TLeg.Heures += Off;
777TLeg.Heures = TLeg.Heures%24;
778if (TLeg.Heures < 0) TLeg.Heures += 24;
779return(TLeg);
780}
781
782
783/* Nouvelle-Fonction */
784HMS TUtoTSolM(HMS TU)
785
786/* Passage TU -> Temps local */
787{
788double t;
789t = HMStoH(TU) + TSolMOff;
790return(HtoHMS(t));
791}
792
793/*
794++
795 Titre Temps Sideral
796| Jour Sideral = 86400 Sec. Sid. ( = 24 H Sid. )
797| Annee Solaire = 365 J 5H 48' 46" =~ 365.2422 J
798| Annee Solaire = 365.2422 Jours Solaire Moyen = 366.2422 Jours Sideraux
799| Jour Sol. Moyem = 24 H Sol. = 24 H Sid. * (366.2422 / 365.2422)
800| Jour Solaire Moyen = 86636.555359 Sec. Sid. = 24 H 3' 56.555359"
801| Heure Solaire Moyen = 3609.85647 Sec. Sid.
802| Sec. Sol. = 1.0027379 Sec. Sid.
803| Le 1er Janvier 1989 0H TU Il etait 6H 42' 30" GMT_Sideral
804--
805*/
806
807/*
808++
809HMS TLegtoTSid(JMA date, HMS TZ)
810 Cette fonction calcule le temps sideral a partir du temps
811 legal et de la date
812GMST_at_0h_UT(JMA date)
813 Cette fonction calcule le temps sideral moyen a greenwich a 0h UT
814--
815*/
816
817/* Nouvelle-Fonction */
818HMS TLegtoTSid(JMA date, HMS TZ)
819
820/* Cette fonction donne le temps sideral a partir du temps legal */
821/* et la date */
822{
823HMS Orgh,TU,TSid;
824JMA Orgd;
825int off;
826long lj;
827double ts,ratio;
828
829/* Le 1er Janvier 1989 0H TU Il etait 6H 42' 30" GMT_Sideral */
830StrgtoJMA("1/1/1989",Orgd);
831StrgtoHMS("6:42:30",Orgh);
832
833
834
835/* Nb de jours / Date Origine */
836lj = JMAtoJ(date) - JMAtoJ(Orgd);
837
838/* On calcule l'heure TU et on tient compte du decalage eventuel de + ou - 1 jour */
839off = TLegOffset(date, TZ);
840TU = TZ;
841TU.Heures -= off;
842if (TU.Heures < 0)
843 {
844 lj -= 1;
845 TU.Heures += 24;
846 }
847if (TU.Heures >= 24)
848 {
849 lj += 1;
850 TU.Heures -= 24;
851 }
852
853
854
855ratio = ((365.*24.) + 5. + (48. / 60.) + (46. / 3600.)) / 24.;
856ratio = (ratio+1) / ratio ;
857
858/* Calcul du decalage en heures du temps Sideral par rapport a */
859/* l'origine choisie : ts */
860
861ts = lj * (24. * (ratio - 1) );
862ts += ( TU.Heures * ratio );
863ts += ( ( TU.Minutes / 60. ) * ratio );
864ts += ( (TU.Secondes / 3600. ) * ratio );
865
866ts += HMStoH(Orgh);
867
868/* On tient compte de l'offset du temps solaire moyen */
869ts += TSolMOff;
870
871TSid = HtoHMS(ts);
872TSid.Heures = TSid.Heures % 24;
873if (TSid.Heures < 0) TSid.Heures += 24;
874
875return(TSid);
876}
877
878
879/* Nouvelle-Fonction */
880HMS GMST_at_0h_UT (JMA date)
881
882/* Cette fonction donne le temps sideral moyen a greenwich a 0h UT */
883/* CMV 5/11/93 */
884{
885HMS Orgh, TSid;
886JMA Orgd;
887long lj, ejs;
888double js, hs, hs0;
889
890/* J H ' " sont les vraies valeurs indiscutables (cf metre etalon etc..) */
891/* Annee Solaire = 365 J 5H 48' 45.2" =~ 365.2421898148 J */
892/* Annee Siderale = 365 J 6H 9' 9.8" =~ 365.2563634259 J */
893/* Le 1er Janvier 1994 0H TU Il sera 6h 41' 40.4305" GMT_Sideral */
894StrgtoJMA("1/1/1994",Orgd);
895StrgtoHMS("6:41:40.4305",Orgh);
896hs0 = HMStoH(Orgh);
897
898/* Nb de jours ecoules depuis la Date Origine */
899lj = JMAtoJ(date) - JMAtoJ(Orgd);
900/* printf("nombre de jours ecoules depuis 1/1/94 %d\n",lj); */
901
902/* Nb de jours sideraux depuis la Date Origine */
903js = (double) lj * 366.2421898148 / 365.2421898148;
904/* Nb d'heures siderales depuis Date Origine (modulo un jour) */
905ejs = (int) js;
906hs = js - (double) ejs;
907hs *= 24.;
908hs += hs0;
909
910TSid = HtoHMS(hs);
911/* printf("nombre d'heures siderales depuis la reference (modulo un jour) %d %d %f\n"
912 ,TSid.Heures,TSid.Minutes,TSid.Secondes); */
913
914return(TSid);
915}
916
917/*
918++
919HMS TUtoTSid(JMA date, HMS TU)
920 Cette fonction calcule le temps sideral GMT a partir du temps
921 universel et de la date GMT.
922int TSidtoTU(JMA date, HMS TS, HMS *TU1, HMS *TU2)
923 Cette fonction calcule les temps universels a partir du temps
924 sideral GMT et de la date GMT. Selon la valeur du TS, il peut y avoir
925 une ou deux possibilites pour le TU (return code).
926--
927*/
928
929/* Nouvelle-Fonction */
930HMS TUtoTSid(JMA date, HMS TU)
931
932/* Cette fonction donne le temps sideral GMT a partir du temps universel */
933/* et la date GMT */
934{
935HMS Orgh,TSid;
936JMA Orgd;
937long lj;
938double ts,ratio;
939
940/* Le 1er Janvier 1989 0H TU Il etait 6H 42' 30" GMT_Sideral */
941StrgtoJMA("1/1/1989",Orgd);
942StrgtoHMS("6:42:30",Orgh);
943
944/* Nb de jours / Date Origine */
945lj = JMAtoJ(date) - JMAtoJ(Orgd);
946
947ratio = ((365.*24.) + 5. + (48. / 60.) + (46. / 3600.)) / 24.;
948ratio = (ratio+1) / ratio ;
949
950/* Calcul du decalage en heures du temps Sideral par rapport a */
951/* l'origine choisie : ts */
952
953ts = lj * (24. * (ratio - 1) );
954ts += ( TU.Heures * ratio );
955ts += ( ( TU.Minutes / 60. ) * ratio );
956ts += ( (TU.Secondes / 3600. ) * ratio );
957
958ts += HMStoH(Orgh);
959
960TSid = HtoHMS(ts);
961TSid.Heures = TSid.Heures % 24;
962if (TSid.Heures < 0) TSid.Heures += 24;
963
964return(TSid);
965}
966
967/* Nouvelle-Fonction */
968int TSidtoTU(JMA date, HMS TS, HMS *TU1, HMS *TU2)
969
970/* Cette fonction donne les temps universels a partir du temps sideral GMT */
971/* et la date GMT */
972{
973int rc = 0;
974HMS Orgh,T0;
975double ratio,delt,deltt,jsecsid;
976char strd[16],strh[16];
977/* char strh1[16]; */
978
979TU1->Heures = TU2->Heures = -1;
980TU1->Minutes = TU2->Minutes = 0;
981TU1->Secondes = TU2->Secondes = 0.;
982
983ratio = ((365.*24.) + 5. + (48. / 60.) + (46. / 3600.)) / 24.;
984ratio = (ratio / (ratio+1));
985jsecsid = 86400./ratio;
986
987/* temps sideral pour la date a 0h TU */
988StrgtoHMS("00:00:00",T0);
989Orgh = TUtoTSid(date,T0);
990/* JMAtoStrg(date,strd); HMStoStrg(T0,strh); HMStoStrg(Orgh,strh1);
991printf("... date=%s tu=%s -> org_ts=%s\n",strd,strh,strh1); */
992
993/* difference de temps sideral / a 0h TU */
994delt = HMStoSec(TS) - HMStoSec(Orgh);
995/* printf("... ts=%f org_ts=%f en sec -> delt=%f (jsecsid=%g)\n"
996 ,HMStoSec(TS),HMStoSec(Orgh),delt,jsecsid); */
997
998if(fabs(delt)<1.e-10) delt=0.;
999if(delt<0.) delt += 86400.;
1000if( fabs(delt-86400.)<1.e-10 ) delt = 86400. - 1.e-10;
1001if(delt<0. || delt>=86400.+1.e-10) {
1002 HMStoStrg(TS,strh);
1003 JMAtoStrg(date,strd);
1004 printf("*** BUG *** TSidtoTU_erreur: delt=%g<0. (d=%s,ts=%s) *******\n"
1005 ,delt,strd,strh);
1006 return(0);
1007}
1008
1009/* difference de temps / a 0h TU 1ere solution */
1010/* printf("... delt=%f en sec sid",delt); */
1011deltt = delt*ratio;
1012*TU1 = SectoHMS(deltt);
1013*TU2 = *TU1;
1014rc++;
1015/* HMStoStrg((*TU1),strh1);
1016printf("... deltt=%f en sec -> TU1=%s\n",deltt,strh1); */
1017
1018/* Y a t-il une seconde solution ? */
1019if(delt>=jsecsid-86400.-1.e-10) return(rc);
1020delt += 86400.;
1021deltt = delt*ratio;
1022*TU2 = SectoHMS(deltt);
1023/* HMStoStrg((*TU2),strh1);
1024printf("... deltt=%f en sec -> TU2=%s\n",deltt,strh1); */
1025rc++;
1026
1027return (rc);
1028}
1029
1030/*
1031++
1032void TSidSetupLaSilla()
1033 Fonction d'initialisation de conversion de temps sideral et legal
1034 pour la ESO-Silla "SetTSolMOff() , SetTLegOff()".
1035| Longitude = 70 deg 43.8 min ouest
1036| Latitude = 29 deg 15.4 min Sud
1037| T Sol. Moyen - TU = -4.715333 Heures
1038| TLegal-TU : -4 Heures de Mars a Octobre, -3 Heures sinon
1039| Changement le 1er dimanche >= 8/03 8/10 a midi (12:00:00)
1040--
1041*/
1042
1043/* Fonction d'initialisation de calcul de temps sideral pour la Silla */
1044void TSidSetupLaSilla()
1045{
1046SetTSolMOff((double)(-4.715333),"La Silla (Chili)");
1047if (SetTLegOff("C -3 87/3 -4 87/10 -3","12:00:00") != 0)
1048 printf(" TSidSetupLaSilla / Erreur SetTLegOff \n");
1049/* Position de l'observatoire */
1050/* SetObsPos("La Silla (Chili)","O 70:43:48","S 29:15:24",2347); */
1051return;
1052}
1053
1054
1055
1056/*
1057++
1058double ToJulianDay(JMA dateTU, HMS hmsTU);
1059 Calcul du jour Julien pour une date TU.
1060 Uniquement valable a partir du 15/10/1582 00:00:00.
1061--
1062*/
1063
1064double ToJulianDay(JMA dateTU, HMS hmsTU)
1065/* Cf Fundamental astronomie, Springer Verlag 2sd ed. cmv 4/12/98 */
1066{
1067double j_dat_0_ut;
1068long int f,g,a;
1069f = (dateTU.Mois>=3) ? dateTU.Annee: dateTU.Annee-1;
1070g = (dateTU.Mois>=3) ? dateTU.Mois: dateTU.Mois+12;
1071a = 2 - (int)((double)f/100.) + (int)((double)f/400.);
1072/* Julian date at 0h UT */
1073j_dat_0_ut = (int)(365.25*(double)f)+(int)(30.6001*(double)(g+1))+dateTU.Jour+a+1720994.5;
1074return j_dat_0_ut + HMStoH(hmsTU)/24.;
1075}
1076
1077/*
1078++
1079int FromJulianDay(double JD,JMA* dateTU, HMS* hmsTU);
1080 Calcul de la date date et l'heure TU a partir d'un jour Julien.
1081 Retourne 0 si succes.
1082--
1083*/
1084
1085int FromJulianDay(double JD,JMA* dateTU, HMS* hmsTU)
1086/* Cf Fundamental astronomie, Springer Verlag 2sd ed. cmv 4/12/98 */
1087{
1088double jd,i,f,a,x,b,c,d,e,h;
1089if(JD<0.) return 1;
1090
1091i = JD - (int)(JD);
1092if(i<0.5) jd = (int)(JD)-0.5; else jd = (int)(JD)+0.5;
1093
1094i = (int)(jd+0.5);
1095f = jd+0.5 - i;
1096if(i<2299161.0) a = i;
1097else {x = (int)((i-1867216.25)/36524.25); a = i+1.+x-(int)(x/4.);}
1098b = a+1524.;
1099c = (int)((b-122.1)/365.25);
1100d = (int)(365.25*c);
1101e = (int)((b-d)/30.6001);
1102
1103x = b-d-(int)(30.6001*e)+f;
1104dateTU->Jour = (int)(x);
1105
1106if(e<13.5) dateTU->Mois = (int)(e-1.); else dateTU->Mois = (int)(e-13.);
1107
1108if(dateTU->Mois>2.5) dateTU->Annee = (int)(c-4716.); else dateTU->Annee = (int)(c-4715.);
1109
1110hmsTU->Heures = hmsTU->Minutes = 0; hmsTU->Secondes = 0.;
1111h = (JD-ToJulianDay(*dateTU,*hmsTU)) * 24.;
1112*hmsTU = HtoHMS(h);
1113
1114return 0;
1115}
1116
1117
1118/*
1119++
1120double StrgtoDegDec(char *strg)
1121 Cette fonction renvoie la valeur decimale en heures
1122 d'un angle specifie sous forme de [-][+]dd:mm:ss.
1123--
1124*/
1125
1126double StrgtoDegDec(char *strg)
1127/* cmv 4/12/98 */
1128{
1129int dd,mm;
1130float ss;
1131char sgn;
1132double deg;
1133
1134dd = mm = 0;
1135ss = 0.;
1136sgn = *strg;
1137if (isdigit(sgn))
1138 sscanf(strg,"%d:%d:%f",&dd,&mm,&ss);
1139else
1140 sscanf(strg+1,"%d:%d:%f",&dd,&mm,&ss);
1141deg = (double)(dd) + (double)(mm)/60. + (double)(ss)/3600.;
1142if (sgn == '-') deg = (-deg);
1143return(deg);
1144}
1145
1146/*
1147++
1148char * DegDectoStrg(double deg, char *strg)
1149 Cette fonction ecrit une valeur en degre decimal sous
1150 forme de [-]dd:mm:ss.
1151--
1152*/
1153
1154char * DegDectoStrg(double deg, char *strg)
1155/* cmv 4/12/98 */
1156{
1157int dd,mm;
1158float ss;
1159char sgn;
1160
1161if (deg < (double)(0.))
1162 { deg = (-deg) ; sgn = '-' ; }
1163else sgn = '+';
1164dd = (int)deg;
1165deg = (deg-(double)(dd)) * (double)60.;
1166mm = (int)deg;
1167ss = (deg-(double)(mm)) * (double)60.;
1168if (sgn == '-')
1169 sprintf(strg,"-%d:%d:%f",dd,mm,ss);
1170else
1171 sprintf(strg,"%d:%d:%f",dd,mm,ss);
1172return(strg);
1173}
1174
1175/*
1176++
1177double EccEarth(JMA dateTU)
1178 Valeur de l'eccentricite de l'orbite terrestre
1179 a la date TU ``dateTU''.
1180--
1181*/
1182double EccEarth(JMA dateTU)
1183/* Cf Fundamental astronomie, Springer Verlag 2sd ed p477 table E10. cmv 9/12/98 */
1184{
1185double JD,T;
1186HMS hmsTU;
1187hmsTU.Heures = 12; hmsTU.Minutes = 0; hmsTU.Secondes = 0.;
1188JD = ToJulianDay(dateTU,hmsTU);
1189/* T en siecles juliens depuis 31/12/1899 a 12h TU --> JD=2415020. */
1190T = (JD-2415020.)/36525;
1191return(0.0167498 - 0.00004258*T - 0.000000137*T*T);
1192}
1193
1194/*
1195++
1196double ObliqEarth(JMA dateTU)
1197 Valeur de l'obliquite de l'orbite terrestre
1198 a la date TU ``dateTU'' (en degres decimaux).
1199--
1200*/
1201double ObliqEarth(JMA dateTU)
1202/* Cf Fundamental astronomie, Springer Verlag 2sd ed. cmv 9/12/98 */
1203/* le 1/1/1900 a 0h00 TU --> jd = 2415020.5 */
1204{
1205double jd;
1206HMS hmsTU;
1207hmsTU.Heures = hmsTU.Minutes = 0; hmsTU.Secondes = 0.;
1208jd = ToJulianDay(dateTU,hmsTU) - 2415020.5;
1209return(23.452294 - 0.00013 *jd/365.25);
1210}
1211
1212/*
1213++
1214double LongEcPerihelie(JMA dateTU)
1215 Retourne la Longitude Ecliptique du perihelie de
1216 l'orbite terrestre a la date TU ``dateTU'' (en degres decimaux).
1217--
1218*/
1219double LongEcPerihelie(JMA dateTU)
1220/* Cf Fundamental astronomie, Springer Verlag 2sd ed p477 table E10. cmv 9/12/98 */
1221{
1222double JD,T;
1223HMS hmsTU;
1224hmsTU.Heures = 12; hmsTU.Minutes = 0; hmsTU.Secondes = 0.;
1225JD = ToJulianDay(dateTU,hmsTU);
1226/* T en siecles juliens depuis 31/12/1899 a 12h TU --> JD=2415020. */
1227T = (JD-2415020.)/36525;
1228/* l = 101d 13' 07.15" + 6171.77" *T + 1.823" *T^2 */
1229return(101.218652778 + 1.71438055556*T + 0.000506388889*T*T);
1230}
1231
1232/*
1233++
1234void EquatToEclip(double a,double d,double* l,double *b,JMA dateTU);
1235 Renvoie les coordonnees ecliptiques ``l,b'' a partir de
1236 coordonnees equatoriales ``a,d'' a la date TU ``dateTU''.
1237 La date permet de corriger la variation annuelle
1238 de l'obliquite terrestre.
1239--
1240*/
1241void EquatToEclip(double a,double d,double* l,double *b,JMA dateTU)
1242/* Cf Fundamental astronomie, Springer Verlag 2sd ed p477 table E10. cmv 9/12/98 */
1243/* l=[0,360[, b=[-90,90] */
1244{
1245double e,sinb,cosb,sinl,cosl;
1246
1247e = ObliqEarth(dateTU)*DegenRad;
1248a *= DegenRad; d *= DegenRad;
1249
1250sinb = sin(d)*cos(e)-cos(d)*sin(e)*sin(a);
1251if(sinb>1.) sinb=1.; if(sinb<-1.) sinb=-1.;
1252*b = asin(sinb);
1253cosb = cos(*b);
1254
1255*l = 0.;
1256if(cosb!=0.) {
1257 sinl = (sin(d)*sin(e)+cos(d)*cos(e)*sin(a))/cosb;
1258 cosl = cos(d)*cos(a)/cosb;
1259 if(cosl>1.) cosl=1.; if(cosl<-1.) cosl=-1.;
1260 *l = acos(cosl);
1261 if(sinl<0.) *l += M_PI;
1262 if(*l>2.*M_PI) *l -= 2.*M_PI;
1263}
1264
1265*b /= DegenRad;
1266*l /= DegenRad;
1267}
Note: See TracBrowser for help on using the repository browser.