#include "dims.h" SUBROUTINE AURHT(ALFA,FLUX,DRIZL,t7,t9,IM,IHEM) ! ! Called by heelis: ! call aurht(t6(3),t8(3),t3(3),t7(3),t9(3),imax,IHEM) ! implicit none C **** CALCULATES HORIZONTAL VARIATION OF AURORAL HEAT SOURCE C **** THIS IS THE NEW VERSION WITH UPDATED ALFA ---- 10/86 ! ! Defines ALFA, FLUX, DRIZL. These are passed back through flowv2, flowv1, ! flowx0, and flowxx to heelis ! #include "params.h" #include "cflowv3.h" #include "ioncr.h" #include "ovalr.h" #include "aurlp.h" #include "phys.h" real :: alfa3,flux3 ! ! Args: integer,intent(in) :: IHEM,im real,intent(out) :: ALFA(IM),FLUX(IM),DRIZL(IM),t7(IM),t9(IM) ! ! Local: integer :: i real :: alfap1,alfap2,pe1,pe2,alfap0,e0p ! ! ALFA0 was defined by tail.f: ! IF (ALFA0 .GT. 0.01) THEN C **** OLD ALFA: ALFA(I) = ALFA0*(1.-RALFA*WK1(I)) DO 1 I=1,IM C **** C **** WK1 = COS(LAMDA) C **** WK1(I)=COS(ATAN2(SIN(ALON(I)-RROTE),COS(ALON(I)-RROTE))) C **** C **** WK2 = AURORAL HALF WIDTH C **** WK2(I)=H0*(1.-RH*COS(ATAN2(SIN(ALON(I)-RROTH), 1 COS(ALON(I)-RROTH))) ) C **** C **** WK3 = COLAT - RAD C **** WK3(I) = COLAT(I) - RRAD(IHEM) C **** OLD ALFA: ALFA(I) = ALFA0*(1.-RALFA*WK1(I)) ALFA(I) = ALFA0*(1.-RALFA*WK1(I)) ! write(6,"(/'aurht: i=',i2,' ihem=',i2)") i,ihem ! write(6,"(' alfa0=',e12.4,' alon(i)=',e12.4)") alfa0,alon(i) ! write(6,"(' rrote=',e12.4,' h0=',e12.4,' rh=',e12.4)") ! | rrote,h0,rh ! write(6,"(' rroth=',e12.4,' rrad(ihem)=',e12.4)") ! | rroth,rrad(ihem) ! write(6,"(' ralfa=',e12.4)") ralfa 1 CONTINUE ! ! ALFA0 <= 0.01: ELSE C **** NEW ALFA (10/86): C **** THE EQUATORWARD DISPLACEMENT OF HARDER ENERGY AROUND C **** 10 MLT IS RD1=RD1A+RD1V*COS(ALON(I)) DO 2 I=1,IM C **** C **** WK1 = COS(LAMDA) C **** WK1(I)=COS(ATAN2(SIN(ALON(I)-RROTE),COS(ALON(I)-RROTE))) C **** C **** WK2 = AURORAL HALF WIDTH C **** WK2(I)=H0*(1.-RH*COS(ATAN2(SIN(ALON(I)-RROTH), 1 COS(ALON(I)-RROTH))) ) C **** C **** WK3 = COLAT - RAD C **** WK3(I) = COLAT(I) - RRAD(IHEM) C **** NEW ALFA (10/86): C **** THE EQUATORWARD DISPLACEMENT OF HARDER ENERGY AROUND C **** 10 MLT IS RD1=RD1A+RD1V*COS(ALON(I)) C ALFA(I,1) = ALFK+ALF1*EXP(-((WK3(I)-(RD1A+RD1V*COS(ALON(I)))) C 1 / RAH1)**2)*EXP(-( (ABS(ALON(I)-RAROT1)-IFIX(ABS(ALON(I)- C 2 RAROT1)/180.)*360.)/RAHT1)**2)+ALF2 * EXP(-( (WK3(I)-RD2) / C 3 RAH2)**2)*EXP(-( (ABS(ALON(I)-RAROT2)-IFIX(ABS(ALON(I)-RAROT2) C 4 /180.)*360.) / RAHT2)**2) C **** C **** ANOTHER NEW ALFA (3/89): C **** ALON IS BETWEEN -180 (0 LT) AND +180. ROT6 IS AROUND C **** -90. AND ROT21 IS AROUND 135. NEED DIFFERENCE 0-180. C **** ALFA(I) = ALFK 1 + ALF6 * EXP(-( (WK3(I)-RD6-RD6V*COS(ALON(I)) )/RH6)**2) * 2 EXP(-( ATAN2(SIN(ALON(I)-RROT6),COS(ALON(I)-RROT6))/RT6 )**2) 3 + ALF21 * EXP(-( (WK3(I)-RD21)/RH21)**2) * 4 EXP(-( ATAN2(SIN(ALON(I)-RROT21),COS(ALON(I)-RROT21))/ 5 RT21)**2) 2 CONTINUE ENDIF ! write(6,"('aurora_heat: lat=',i3,' alfa=',/,(6e12.4))") ! | j,alfa ! DO I=1,IM FLUX(I)=E0*(1.-REE*WK1(I))*EXP(-( WK3(I) /WK2(I))**2) | / (2. * ALFA(I) * 1.602E-9) DRIZL(I)=EXP(-(( WK3(I) +ABS( WK3(I) ))/(2.*H0))**2) t7(i) = ALFA20*(1.-RALFA2*WK1(I)) t9(i) = E20*(1.-RE2*WK1(I))*EXP(-( WK3(I) /WK2(I))**2) | / (2. * t7(i) * 1.602E-9) c c Low-energy proton inputs: ALFAP1 = 5. ALFAP2 = 15. C PE1 = 0.15 C PE2 = 0.4 PE1 = 1.E-20 PE2 = 1.E-20 ALFAP0 = 0.5*(ALFAP1+ALFAP2) alfa_lp(i+2) = ALFAP0*(1-RALFA2*WK1(I)) E0P = 0.5*(PE1+PE2) eflux_lp(i+2) = E0P*(1.-REE*WK1(I))*exp(-(wk3(i)/wk2(i))**2) qteaur(i+2) = -7.e+8*exp(-(wk3(i)/wk2(i))**2) enddo ! write(6,"('aurora_heat: lat=',i3,' flux=' ,/,(6e12.4))") j,flux ! write(6,"('aurora_heat: lat=',i3,' drizl=',/,(6e12.4))") j,drizl ! write(6,"('aurora_heat: lat=',i3,' alfa2=',/,(6e12.4))") j,t7 ! write(6,"('aurora_heat: lat=',i3,' flux2=',/,(6e12.4))") j,t9 c c Periodic points for alfa_lp and eflux_lp: alfa_lp(1) = alfa_lp(im+1) alfa_lp(2) = alfa_lp(im+2) alfa_lp(im+3) = alfa_lp(3) alfa_lp(im+4) = alfa_lp(4) eflux_lp(1) = eflux_lp(im+1) eflux_lp(2) = eflux_lp(im+2) eflux_lp(im+3) = eflux_lp(3) eflux_lp(im+4) = eflux_lp(4) qteaur(1) = qteaur(im+1) qteaur(2) = qteaur(im+2) qteaur(im+3) = qteaur(3) qteaur(im+4) = qteaur(4) RETURN END C