solter15.f
4.63 KB
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
subroutine solter15 (iyear,imonth,iday,ihour,imin,isec,year,
> alfag,alfas,deltas,obliq,tetdip,phidip,ier)
c*
c***********************************************************************
c*
c* "Copyright [c] CNES 98 - tous droits reserves"
c* **********************************************
c*
c*PRO MAGLIB
c*
c*VER 2017.02.24 - V 1.0
c*
c*AUT spec. CNES - JC KOSIK - octobre 2001
c*AUT port. CISI
c*AUT adapt. CISI
c*
c*ROL Theme : Astronomie et calculs d'orbite
c*ROL Calcul de l'ascension droite de Greenwich.
c*ROL Calcul de l'ascension droite et de la declinaison
c*ROL du soleil.
c*ROL Calcul de la position geographique du point ou le dipole
c*ROL coupe l'hemisphere nord.
c*
c*PAR iyear (I) : annee de la date (annee >= 2015)
c*PAR imonth (I) : mois de la date
c*PAR iday (I) : jour de la date
c*PAR ihour (I) : heures de la date
c*PAR imin (I) : minutes de la date
c*PAR isec (I) : secondes de la date
c*
c*PAR year (O) : annee fractionnaire (>= 2015)
c*
c*PAR alfag (O) : ascension droite de Greenwich (radians)
c*
c*PAR alfas (O) : ascension droite du soleil (radians)
c*PAR deltas (O) : declinaison du soleil (radians)
c*
c*PAR obliq (O) : obliquite de l'ecliptique (radians)
c*
c*PAR tetdip (O) : colatitude geocentrique du point ou le dipole
c*PAR : coupe l'hemisphere nord (radians)
c*PAR phidip (O) : longitude geocentrique du point ou le dipole
c*PAR : coupe l'hemisphere nord (radians)
c*
c*PAR ier (O) : code de retour
c*
c*NOT ier : sans objet
c*
c*INF utilise : julg, sun
c*
c*HST version 1.0 - 2017.02.24 - Utilisation du champ IGRF15 (annee fractionnaire
c*HST entre 2015 et 2020)
c*
c***********************************************************************
c*
implicit none
c
c ---------------------------------
c*FON Declaration identificateur rcs_id
c ---------------------------------
c
character rcs_id*100
c
c --------------------------
c*FON Declaration des parametres
c --------------------------
c
integer iyear, imonth, iday, ihour, imin, isec
double precision year
double precision alfag, alfas, deltas
double precision obliq, tetdip, phidip
integer ier
c
c ----------------------------------
c*FON Declaration des fonctions externes
c ----------------------------------
c
external julg
integer julg
c
c ---------------------------------
c*FON Declaration des variables locales
c ---------------------------------
c
integer julc,jul1,nujour
c
double precision fjour,dt,g10,g11,c11,h11,gdt
c*LOC Variables de travail intermediaires
c
SAVE
c
c ---------------------------------
c*FON Affectation identificateur rcs_id
c ---------------------------------
c
data rcs_id /"
>$Id$"/
c
c ******************
c Debut de programme
c ******************
c
ier = 0
c
c -----------------------------
c*FON Calcul de la date en secondes
c -----------------------------
c
fjour = dble((ihour * 60 + imin) * 60 + isec) / 86400.d0
c
c --------------------------------
c*FON Calcul du jour julien CNES 1950
c --------------------------------
c
julc = julg(iday,imonth,iyear,ier)
c
jul1 = julg(1,1,iyear,ier)
c
c ----------------------------------------------------------
c*FON Calcul du numero du jour compte a partir du 1er janvier de
c*FON l'annee courante
c ----------------------------------------------------------
c
nujour = julc - jul1 + 1
c
c -----------------------------------------------------------
c*FON Calcul de la position du soleil dans le repere equatorial a
c*FON gamma 50
c -----------------------------------------------------------
c
call sun(iyear,nujour,fjour,alfag,alfas,deltas,obliq,ier)
c
year = dble(iyear) + dble(nujour-1) / 365.d0
c
c -------------------------------------------------------------
c*FON Calculs intermediaires pour obtenir les coordonnees du dipole
c -------------------------------------------------------------
c
dt = year - 2010.d0
g10 = +29442.0d0 - 10.3d0 * dt
g11 = +1501.0d0 - 18.1d0 * dt
h11 = -4997.1d0 + 26.6d0 * dt
c
c11 = sqrt(g11**2 + h11**2)
gdt = sqrt(g11**2 + h11**2 + g10**2)
c
c ---------------------------------------------------------------
c*FON Calcul de la longitude et de la colatitude du dipole en radians
c ---------------------------------------------------------------
c
phidip = atan(h11 / g11)
tetdip = asin(c11 / gdt)
c
c ****************
c Fin de programme
c ****************
c
return
end