k899k.f
3.79 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
subroutine k899k (tilt,facrc,c20,rre,thet,phi,bre,bte,bpe,ier)
c*
c***********************************************************************
c*
c* "Copyright [c] CNES 98 - tous droits reserves"
c* **********************************************
c*
c*PRO MAGLIB
c*
c*VER 99.03.31 - V 1.0
c*VER 01.06.05 - V 2.0
c*VER 03.01.06 - V 2.1
c*
c*AUT spec. CNES - JC KOSIK - janvier 1998
c*AUT port. CISI
c*
c*ROL Theme : Modeles de champs magnetiques
c*ROL Calcul du courant de l'anneau du champ magnetique.
c*
c*PAR tilt (I) : angle de tilt (radians)
c*
c*PAR facrc (I) : coefficient du courant de l'anneau
c*PAR c20 (I) : coefficient du courant de l'anneau
c*
c*PAR rre (I) : distance radiale (kilometres ou rayons terrestres)
c*PAR thet (I) : colatitude geocentrique GSM (radians)
c*PAR phi (I) : colatitude geocentrique GSM (radians)
c*
c*PAR bre (O) : distante radiale geocentrique (idem input)
c*PAR bte (O) : colatitude geocentrique (radians)
c*PAR bpe (O) : colatitude geocentrique (radians)
c*
c*PAR ier (O) : code de retour
c*
c*NOT ier : sans objet
c*
c*INF utilise : sans objet
c*
c*HST version 1.0 - 99.03.31 - creation de la maglib au CDPP
c*HST version 2.0 - 01.06.05 - correction de commentaires de code
c*HST version 2.1 - 03.01.06 - corrections en compilation avec g77
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
double precision tilt, facrc, c20, rre, thet, phi
double precision bre, bte, bpe
integer ier
c
c ---------------------------------
c*FON Declaration des variables locales
c ---------------------------------
c
double precision br,bt,bp
c*LOC br,bt,bp : composantes radiale, tangentielle et azimuthale du champ
c
double precision brtilt,bttilt,bptilt
c*LOC brtilt,bttilt,bptilt : composantes radiale, tangentielle et azimuthale
c
double precision c10,c21,xk1,xk2,xk2t,r2,r3
double precision fons1,fons2,fons2t,ct,st,s2t,c2t,cp,sp
double precision dfons1, dfons2, dfons2t
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
c10 = -1.5d0
c21 = 0.11d0
c
xk1 = 0.04d0
xk2 = 0.01d0
xk2t = 0.005d0
c
r2 = rre * rre
r3 = r2 * rre
c
fons1 = r3 * dexp(-xk1 * r2)
fons2 = r3 * dexp(-xk2 * r2)
fons2t = r2 * dexp(-xk2t * r2)
c
ct = dcos(thet)
st = dsin(thet)
s2t = 2.d0 * st * ct
c2t = ct * ct - st * st
cp = dcos(phi)
sp = dsin(phi)
c
dfons1 = (4.d0 - 2.d0 * xk1 * r2) * fons1 / rre
dfons2 = (4.d0 - 2.d0 * xk2 * r2) * fons2 / rre
dfons2t = (3.d0 - 2.d0 * xk2t * r2) * fons2t / rre
c
br = 2.d0 * c10 * fons1*ct
> / rre + 6.d0 * c21 * ct * st * cp * fons2 / rre
bt = -c10 * dfons1 * st + 1.732d0 * c21 * c2t * cp * dfons2
bp = -1.732d0 * c21 * ct * sp * dfons2
c
br = br * facrc
bt = bt * facrc
bp = bp * facrc
c
brtilt = +c20 * (6.d0 * c2t - 2.d0) * dsin(tilt)
> * facrc * fons2t / rre
bttilt = -c20 * s2t * dsin(tilt) * dfons2t * facrc
bptilt = 0.d0
c
bre = br * dcos(tilt) + brtilt
bte = bt * dcos(tilt) + bttilt
bpe = bp * dcos(tilt) + bptilt
c
c ****************
c Fin de programme
c ****************
c
999 continue
return
end