subroutine kk97kp (indval,tilt,xgsm,ygsm,zgsm, > bxt,byt,bzt,bt,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.05.30 - 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 champ externe dans un systeme en coordonnees c*ROL solaires magnetospheriques. Les unites sont les radians, c*ROL rayons terrestres et nanoteslas. c c*PAR indval (I) : indice geomagnetique : niveau d'amplitude du champ c* c*PAR tilt (I) : angle de tilt (radians) c* c*PAR xgsm (I) : coordonnee solaire magnetospherique x (rayons terrestres) c*PAR ygsm (I) : coordonnee solaire magnetospherique y (rayons terrestres) c*PAR zgsm (I) : coordonnee solaire magnetospherique z (rayons terrestres) c* c*PAR bxt (O) : composante solaire magnetospherique x du champ c*PAR : externe (nanoteslas) c*PAR byt (O) : composante solaire magnetospherique y du champ c*PAR : externee (nanoteslas) c*PAR bzt (O) : composant solaire magnetospherique z du champ c*PAR : externe (nanoteslas) c* c*PAR bt (O) : totalite du champ (nanoteslas) c* c*PAR ier (O) : code de retour c* c*NOT indval : 1 ---> Kp = 0 , 0+ c*NOT indval : 2 ---> Kp = 1- , 1 , 1+ c*NOT indval : 3 ---> Kp = 2- , 2 , 2+ c*NOT indval : 4 ---> Kp = 3- , 3 , 3+ c*NOT indval : 5 ---> Kp = 4- , 4 , 4+ c*NOT indval : 6 ---> Kp > 5- c* c*NOT ier : sans objet c* c*INF utilise : carsp, k899k, vspvcar, tailts1k, bdbspst2k c*INF utilise : cylind2k c* c*HST version 1.0 - 99.03.31 - creation de la maglib au CDPP c*HST version 2.0 - 01.05.30 - 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,xgsm,ygsm,zgsm double precision bxt,byt,bzt,bt integer indval,ier c c --------------------------------- c*FON Declaration des variables locales c --------------------------------- c integer i,j,m c*LOC i,j,m : indices de boucles c double precision gg(6,6) c*LOC gg : coefficients du champ lointain c double precision rb c*LOC rb : distance radiale geocentrique c double precision xn,bn,delb,s c*LOC xn,bn,delb,s : coefficients du modele de la queue c double precision facrc,c20 c*LOC facrc,c20 : coefficients du courant de l'anneau c double precision tx(5,17) c*LOC tx : tx(1) ..to tx(9) gg(i,j) coefficients, tx(10) = rb c*LOC tx(11).....tx(14) coefficients xn, bn, delb, s c*LOC tx(15),tx(16) facrc, c20 c*LOC tx(17) = bb c double precision bxts,byts,bzts c*LOC bxts,byts,bzts : composantes solaires magnetospheriques en x, y et z c double precision bb c*LOC bb : coefficient du systeme de courants en retour c double precision bmxe,bmye,bmze c*LOC bmxe,bmye,bmze : composantes xgsm, ygsm et zgsm du champ lointain c double precision bcylx,bcyly,bcylz c*LOC bcylx,bcyly,bcylz : composantes du champ magnetique induit par les c*LOC : courants en retour en x, y et z c double precision xsm,ysm,zsm c*LOC xsm,ysm,zsm : coordonnees solaires magnetiques cartesiennes c double precision rsm,thetsm,phism c*LOC rsm,thetsm,phism : composantes solaires magnetiques spheriques c double precision brsm,btsm,bpsm c*LOC brsm,btsm,bpsm : vecteur du courant de l'anneau c double precision bxsm,bysm,bzsm c*LOC bxsm,bysm,bzsm : coordonnees solaires magnetiques cartesiennes c*LOC : du courant de l'anneau c double precision bxgsm,bygsm,bzgsm c*LOC bxgsm,bygsm,bzgsm : coordonnees GSM cartesiennes c*LOC : du courant de l'anneau c double precision rgsm,thgsm,phgsm c*LOC rsm,thetsm,phism : composantes GSM spheriques c double precision xmax,hmax,xgsmo,rcyl,xk,delz,zgsk c*LOC xmax,hmax,xgsmo : variables intermediaires pour le calcul c*LOC rcyl,xk,delz,zgsk : du champ magnetique c SAVE c c --------------------------------- c*FON Affectation identificateur rcs_id c --------------------------------- c data rcs_id /" >$Id$"/ c c -------------------------- c*FON Affectation des constantes c -------------------------- c data (tx(1,m),m=1,17)/ > -10.98617d0, -8.50034d0, 1.34113d0, -2.48786d0, 0.30964d0, > -0.42732d0, -0.02118d0, 0.01296d0, -0.03237d0, 12.d0, > -7.d0, 50.d0, 45.d0, 70.d0, 1.d0, > -0.6d0, 25.d0/ data (tx(2,m),m=1,17)/ > -9.91009d0, -9.34152d0, 1.68348d0, -2.85354d0, 0.43415d0, > -0.47281d0, -0.02807d0, 0.02365d0, -0.03232d0, 12.d0, > -7.d0, 40.d0, 30.d0, 70.d0, 1.3d0, > -0.6d0, 30.d0/ data (tx(3,m),m=1,17)/ > -12.37934d0, -10.70398d0, 1.86899d0, -3.45258d0, 0.52744d0, > -0.59398d0, -0.02654d0, 0.03671d0, -0.04219d0, 12.d0, > -7.d0, 40.d0, 30.d0, 70.d0, 1.5d0, > -0.4d0, 35.d0/ data (tx(4,m),m=1,17)/ > -10.32178d0, -35.92150d0, 25.41001d0, -38.68030d0, 22.99485d0, > -20.01264d0, -3.43988d0, 5.45664d0, -3.66522d0, 36.d0, > -7.d0, 40.d0, 30.d0, 70.d0, 1.8d0, > -0.4d0, 50.d0/ data (tx(5,m),m=1,17)/ > -3.03703d0, 0.35692d0, 0.12408d0, 1.13281d0, 0.10293d0, > 0.35461d0, 0.00433d0, 0.01610d0, 0.03338d0, 10.d0, > -7.d0, 50.d0, 40.d0, 70.d0, 2.d0, > 0.2d0, 60.d0/ c c ****************** c Debut de programme c ****************** c ier = 0 c c ------------------------------------------------ c*FON Les coefficients gg(i,j) sont initialises a zero c ------------------------------------------------ c do 10 i = 2,6 do 20 j = 1,i gg(i,j) = 0.d0 20 continue 10 continue c c ------------------------------ c*FON Coefficients du champ lointain c ------------------------------ c gg(2,1) = tx(indval,1) gg(2,2) = 0.d0 gg(3,1) = 0.d0 gg(3,2) = tx(indval,2) gg(3,3) = 0.d0 gg(4,1) = tx(indval,3) gg(4,2) = 0.d0 gg(4,3) = tx(indval,4) gg(4,4) = 0.d0 gg(5,1) = 0.d0 gg(5,2) = tx(indval,5) gg(5,3) = 0.d0 gg(5,4) = tx(indval,6) gg(5,5) = 0.d0 gg(6,1) = tx(indval,7) gg(6,2) = 0.d0 gg(6,3) = tx(indval,8) gg(6,4) = 0.d0 gg(6,5) = tx(indval,9) gg(6,6) = 0.d0 c rb = tx(indval,10) c c ---------------------------------- c*FON Coefficients du modele de la queue c ---------------------------------- c xn = tx(indval,11) bn = tx(indval,12) delb = tx(indval,13) s = tx(indval,14) c c ------------------------------------ c*FON Coefficients ddu courant de l'anneau c ------------------------------------ c facrc = tx(indval,15) c20 = tx(indval,16) c c -------------------------------------------- c*FON Coefficient du systeme de courants en retour c -------------------------------------------- c bb = tx(indval,17) c c ------------------------------------------- c*FON Tous les composants sont initialises a zero c ------------------------------------------- c bxts = 0.d0 byts = 0.d0 bzts = 0.d0 c bmxe = 0.d0 bmye = 0.d0 bmze = 0.d0 c bcylx = 0.d0 bcyly = 0.d0 bcylz = 0.d0 c c ----------------------------------------- c*FON Calcul du champ k89 (courant de anneau) : c ----------------------------------------- c c ---------------------------------------------------------- c*FON 1) Les coordonnees GSM sont transformees en coordonnees SM c ---------------------------------------------------------- c xsm = xgsm * dcos(tilt) - zgsm * dsin(tilt) ysm = ygsm zsm = xgsm * dsin(tilt) + zgsm * dcos(tilt) c c ------------------------------------------------------- c*FON 2) Les coordonnees cartesiennes SM sont transformees en c*FON coordonnees spheriques SM c ------------------------------------------------------ c call carsp(xsm,ysm,zsm,rsm,thetsm,phism,ier) c c --------------------------------------------------- c*FON 3) Le champ courant de l'anneau est calcule dans le c*FON systeme de coordonnees speriques SM c --------------------------------------------------- c call k899k (tilt,facrc,c20,rsm,thetsm,phism,brsm,btsm,bpsm,ier) c c ------------------------------------------------- c*FON 4) Les composants spheriques SM sont transformees c*FON en composants cartesiens SM c ------------------------------------------------- c call vspvcar(thetsm,phism,brsm,btsm,bpsm,bxsm,bysm,bzsm,ier) c c ----------------------------------------------------- c*FON 5) Cette contribution est transformee dans le systeme c*FON de coordonnees GSM c ----------------------------------------------------- c bxgsm = bxsm * dcos(tilt) + bzsm * dsin(tilt) bygsm = bysm bzgsm = -bxsm * dsin(tilt) + bzsm * dcos(tilt) c c --------------------------------------------------------------- c*FON Le champ queue est calcule en coordonnees systeme GSM. C'est le c*FON modele de queue 1982 Tsyganenko avec des parametres modifies c --------------------------------------------------------------- c call tailts1k(tilt,xn,bn,delb,s,xgsm,ygsm,zgsm,bxts,byts,bzts,ier) c c ------------------------------------- c*FON Calcul des coordonnees spheriques GSM c ------------------------------------- c call carsp(xgsm,ygsm,zgsm,rgsm,thgsm,phgsm,ier) c c ------------------------------------------ c*FON Calcul du champ magnetique lointain en GSM c ------------------------------------------ c call bdbspst2k(tilt,rb,gg,rgsm,thgsm,phgsm,bmxe,bmye,bmze,ier) c c ----------------------------------------------------------------- c*FON Calcul du champ magnetique a partir des courants en retour en GSM c ----------------------------------------------------------------- c xmax = xn * cos(tilt) hmax = dabs(xn) * sin(tilt) xgsmo = 0.d0 rcyl = 30.d0 xk = 0.05d0 c if (xgsm .ge. xmax) then c delz = xgsm * hmax / xmax zgsk = zgsm - delz call cylind2k(xk,bb,xgsmo,rcyl,xgsm,ygsm,zgsk, > bcylx,bcyly,bcylz,ier) c else if (xgsm .lt. xmax) then c delz = hmax zgsk = zgsm - delz call cylind2k(xk,bb,xgsmo,rcyl,xgsm,ygsm,zgsk, > bcylx,bcyly,bcylz,ier) c endif c c ----------------------------------------------------------------- c*FON Somme des contributions au champ magnetique a partir de l'anneau, c*FON de la queue, du champ lointain et des courants en retour c ----------------------------------------------------------------- c bxt = bxgsm + bxts + bmxe + bcylx byt = bygsm + byts + bmye + bcyly bzt = bzgsm + bzts + bmze + bcylz c c ------------------------------- c*FON Calcul du module du champ total c ------------------------------- c bt = sqrt(bxt * bxt + byt * byt + bzt * bzt) c c **************** c Fin de programme c **************** c return end