c------------------------------------------------------ c to test different SW parameterizations c------------------------------------------------------ open (11,file='sw.res',status='unknown') c>>>>>>> definitions of atmospheric variables <<<<<<<<<< c input: c TA - air temperature, C real*4 c ez - water vapor pressure, mb real*4 c cn - total fractinal cloud cover (%/10) real*4 c cl - low-level fractional cloud cover real*4 c cn8 - total octa cloud cover real*4 C IE - =0, if two levels cloudiness INTEGER C =1, if only total cloudiness c output: c Qs - incoming SW radiation, W/m*m real*4 c Rpogl - net SW (Malevsky), W/m*m real*4 c Rdob - net SW (Dobson), W/m*m real*4 c c>>>>>>> definition of astromomical variables <<<<<<<<<< c input: c ifebr - the number of days in February c in the current year integer*4 c rmo - month real*4 c rlat - latitude real*4 c rlon - longtitude real*4 c da - date real*4 c rhr - hour, UTC real*4 c output: c h - solar altitude real*4 c albedo - albedo (Payne, 1971) real*4 c------------------------------------------------------ IFEWR=28 rlat=23. rlon=12. rmo=6. da=28. hr=12.6 ta=15. ez=15. cn=9.5 cl=1. cn8=7.6 c----------- SW radiation ---------------------- print *, '1' c>>>>>>>>> Dobson - instantaneous <<<<<<<<<<<<<<<<<<<<<< call albpy (rmo,rlat,albedo) print *,albedo call SALT1 (da,rmo,rlat,hr,IFEWR,H) print *, 'h', h call RSWd8 (H,cn8,albedo,rdob) print *, '2' c>>>>>>>> Malevsky - instantaneous<<<<<<<<<<<<<<<<<<<<<<<<<<<<< call SALT1 (da,rmo,rlat,hr,IFEWR,H) call RSW (H,TA,EZ,cl,cn,QS,RPOGL,1) print *, '3' c---------------------------------------------- write (11,1000) H, albedo write (11,1001) rdob write (11,1002) rpogl 1000 format ('solar alt. albedo',2f12.3) 1001 format ('Dobson',f12.3) 1002 format ('Malevsky',f12.3) stop end C******************************************************************** C s h o r t w a v e r a d i a t i o n C monthly mean values C Reed - 1977 C******************************************************************** C AM - month REAL " C SR - latitude REAL " C ANO- common cloudness percent/10 REAL " C RPOGL-accumulate radiation W/M**2, REAL " C******************************************************************* SUBROUTINE RSWisi (AM,SH,ANO,RPOGL) DIMENSION Amon(12),Ash(9),ALB(9,12) DATA Amon/1.,2.,3.,4.,5.,6.,7.,8.,9.,10.,11.,12./ DATA Ash/-80.,-70.,-60.,-50.,-40.,-30.,-20.,-10.,0./ DATA ALB/0.41,0.41,0.28,0.11,0.10,0.09,0.07,0.07,0.06, * 0.41,0.41,0.12,0.10,0.09,0.07,0.06,0.06,0.06, * 0.33,0.15,0.09,0.08,0.07,0.06,0.06,0.06,0.06, * 0.14,0.10,0.07,0.07,0.07,0.06,0.06,0.06,0.06, * 0.10,0.08,0.07,0.06,0.06,0.06,0.06,0.06,0.06, * 0.09,0.07,0.07,0.06,0.06,0.06,0.06,0.06,0.06, * 0.08,0.07,0.06,0.06,0.06,0.06,0.06,0.06,0.06, * 0.08,0.09,0.07,0.07,0.06,0.06,0.06,0.06,0.06, * 0.12,0.10,0.07,0.07,0.07,0.06,0.06,0.06,0.06, * 0.25,0.25,0.10,0.08,0.08,0.07,0.06,0.06,0.06, * 0.25,0.25,0.16,0.11,0.10,0.08,0.07,0.06,0.06, * 0.44,0.44,0.44,0.12,0.11,0.09,0.07,0.07,0.06/ ano1=ano/1.25 c********************************************************************* c albedo calculation -- albedo = F (solar high, total cloud cover) c********************************************************************* sh1=sh*(-1.) CALL INTERP(AM,SH1,Ash,Amon,9,12,ALB,ALP) c********************************************************************* c incoming shortwave radiation c********************************************************************* id=15 im=int(am) call Juld90(im,id,itag) RPOGL=QREED(itag,ano1,sh,alp) RETURN END SUBROUTINE albpy (AM,SH,albedo) DIMENSION Amon(12),Ash(9),ALB(9,12) DATA Amon/1.,2.,3.,4.,5.,6.,7.,8.,9.,10.,11.,12./ DATA Ash/-80.,-70.,-60.,-50.,-40.,-30.,-20.,-10.,0./ DATA ALB/0.41,0.41,0.28,0.11,0.10,0.09,0.07,0.07,0.06, * 0.41,0.41,0.12,0.10,0.09,0.07,0.06,0.06,0.06, * 0.33,0.15,0.09,0.08,0.07,0.06,0.06,0.06,0.06, * 0.14,0.10,0.07,0.07,0.07,0.06,0.06,0.06,0.06, * 0.10,0.08,0.07,0.06,0.06,0.06,0.06,0.06,0.06, * 0.09,0.07,0.07,0.06,0.06,0.06,0.06,0.06,0.06, * 0.08,0.07,0.06,0.06,0.06,0.06,0.06,0.06,0.06, * 0.08,0.09,0.07,0.07,0.06,0.06,0.06,0.06,0.06, * 0.12,0.10,0.07,0.07,0.07,0.06,0.06,0.06,0.06, * 0.25,0.25,0.10,0.08,0.08,0.07,0.06,0.06,0.06, * 0.25,0.25,0.16,0.11,0.10,0.08,0.07,0.06,0.06, * 0.44,0.44,0.44,0.12,0.11,0.09,0.07,0.07,0.06/ c********************************************************************* c albedo calculation -- albedo = F (solar high, total cloud cover) c********************************************************************* if (sh.gt.0.) sh1=sh*(-1.) CALL INTERP(AM,SH1,Ash,Amon,9,12,ALB,albedo) RETURN END cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc c P R A D . F O R cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc C Some English comments: c c This is a part of a package with radiation subroutines c written by hj isemer. c c It essentially contains solar radiation subroutines c for the methods of c REED, 1977, see Journal of Phys. Oceanogr. c and c DOBSON and SMITH, 1988 see Quart.Journ.Royal Met.Soc. c c and additional auxiliary routines c--------------------------------------------------------------------- C 1. FUNCTION QREED BERECHNET DAS TAGESMITTEL C DAZU MUSS DAS MITTEL DER BEWOELKUNG BERECHNET WERDEN, ABER C NUR FUER DIE STUNDEN AM TAG. DIES KANN MIT DER ROUTINE C CLOUDA GEMACHT WERDEN. C 2. FUNCTION QDOBSON BERECHNET DIE STRAHLUNG FUER STUENDLICHE C WOLKENBEOBACHTUNGEN NACH DEM LINEAREN OKTAMODELL VON C DOBSON UND SMITH. C 3. SUBROUTINE QDODAY BERECHNET DAS TAGESMITTEL NACH DOBSON UND C SMITH, DIES KANN MIT DEM ERGEBNIS VON QREED VERGLICHEN WERDEN. C C QDODAY UND CLOUDA MUESSTEN NOCH GETESTET WERDEN. C C DIE ALBEDO DER OZEANOBERFLaeCHE KANN NACH PAYNE, 1972 ALS C KLIMATOLOGISCHER MITTELWERT (ETWA 0.07) ODER WENN MoegLICH C AUS ANDEREN MESSUNGEN EINGESETZT WERDEN. C ZUR BERECHNUNG DER EINFALLENDEN KURZWELLIGEN STRAHLUNG NULL C SETZEN! C C DIE UEBRIGEN ROUTINEN WERDEN Z.T. AUFGERUFEN. C C PRAD.FOR IST UNTER MS FORTRAN ALS PRAD.LIB ALS C BINDEFAEHIGE BIBLIOTHEK VORHANDEN. C C HJ ISEMER, SEPTEMBER 1989 c****************************************************************** FUNCTION SI0(JD) C **************** C Calculates the 'solar constant' as a function of the c Julian day of the year c input : JD is the Julian day c INTEGER JD PI=3.1415927 SOLAR=1368 S1=COS(PI*2*(JD-21)/365)*0.035 SI0=SOLAR*(1+S1) END FUNCTION SUNZEN(PHI,LAM,MON,TAG,ZGMT) C ************************************* c c calculates the zenith angle of the sun as a function of c latitude: PHI c longitude: LAM c month : MON c day in the month : TAG c time of the day (GMT!) : ZGMT REAL PHI,LAM,ZGMT INTEGER MON,TAG DELTA=RAD(DEKLIN(MON,TAG)) H=STUWI(OZEIT(LAM,ZGMT)) XLA=RAD(PHI) S1=(SIN(XLA)*SIN(DELTA))+(COS(XLA)*COS(DELTA)*COS(H)) S2=ACOS(S1) SUNZEN=WI(S2) END C FUNCTION DEKLIN(MON,TAG) C ************************ c calculates the deklination angle of the sun as a function of c the month : MON c the day in the month : TAG INTEGER MON,TAG PI=3.1415927 CALL JULD90(MON,TAG,JULD) SDEKLIN=SIN(RAD(23.45))*SIN((2*PI*(JULD-80)/365)) DEKLIN=WI(ASIN(SDEKLIN)) END C FUNCTION DEKLIJ(JulD) C ************************ c calculates the deklination angle of the sun as a function of c the Julian Day : JulD INTEGER JD PI=3.1415927 SDEKLIN=SIN(RAD(23.45))*SIN((2*PI*(JulD-80)/365)) DEKLIJ=WI(ASIN(SDEKLIN)) END C REAL FUNCTION OZEIT(LAMDA,GMTZEIT) C ************************************ REAL GMTZEIT,LAMDA C BERECHNUNG DER LOKALEN UHRZEIT AUS GEOGRAPHISCHER LAENGE UND C GMT UHRZEIT , LAMDA NEGATIV FUER OESTLICHE LAENGE c c calculates the local time (in full hours only) as a function of c longitude : LAMDA c GMT time : GMTZEIT OZEIT=GMTZEIT+(LAMDA/360*24) END C FUNCTION STUWI(ZEIT) C ******************** REAL ZEIT PI=3.1415927 STUWI=(ZEIT-12)/24*2*PI END C FUNCTION QREED(ITAG,CNTOT,RPHI,ALBEDO) C ******************************* INTEGER ITAG REAL CNTOT,RPHI c calculates daily solar radiation after REED (1977), J.Phys.Oc.7,482 c as a function of c Julian day : ITAG c mean daily total cloud cover : CNTOT (in octas: CNTOT=5 means 5/8) c latitude : RPHI c albedo of the sea surface : ALBEDO C BERECHNUNG VON TAGESMITTELWERTEN (!) C DER VOM OZEAN ABSORBIERTEN C KURZWELLIGEN EINSTRAHLUNG BEI BEWOELKTEM HIMMEL NACH C REED, 1977 QREED IN (W/M**2) C C EINGABE: C ITAG : JULIANISCHER TAG IM JAHR C CNTOT : MITTLERE GESAMTBEDECKUNG DURCH WOLKEN AN DIESEM TAG C IN 1/8, ZB CNTOT=5.5 BEDEUTET 5.5/8. BEDECKUNG C NUR AUS DEN BEOBACHTUNGEN BERECHNEN, AN DENEN DIE C SONNE UEBER DEM HORIZONT STEHT C SIEHE ROUTINE CLOUDA WEITER UNTEN ! C RPHI : GEOGRAPHISCHE BREITE C ALBEDO: ALBEDO DER OZEANOBERFLAECHE CALL SUNALT(ITAG,RPHI,SUNMID) print *, sunmid CALL Q0SEBE(RPHI,ITAG,QCLEAR) QREED=QKCLOU(CNTOT,SUNMID,QCLEAR)*(1.-ALBEDO) RETURN END C FUNCTION QDOBSON(PHI,JD,RSTD,ICL,ALBEDO) C ********************************* c hourly solar radiation at the sea surface from the linear okta modell c of DOBSON and SMITH, see QJRMS, 1988 c as a function of c latitude : PHI c Julian day: JD c hour of the day (0-23): RSTD c total cloud cover in oktas (0-8): ICL c albedo of the sea surface : ALBEDO C C BERECHNUNG DER VOM OZEAN ABSORBIERTEN C KURZWELLIGEN EINSTRAHLUNG NACH C DEM LINEAREN OKTA MODELL VON DOBSON & SMITH, QJRMS 1988 C AUS DATEN DES OWS PAPA IM PAZIFIK C C EINGABE: C C PHI : GEOGR.BREITE C JD : JULIANISCHER TAG IM JAHR C RSTD : STUNDE AM TAG (0 BIS 23) C ICL : GESAMTBEDECKUNG IN ACHTELN, DARF NUR 0 BIS 8 SEIN C ALBEDO: ALBEDO DER OZEANOBEFLAECHE REAL A(0:9),B(0:9) C KOEFFIZIENTEN DES OKTA MODELLS DATA A/.4,.517,.474,.421,.38,.35,.304,.23,.106,.134/ DATA B/.386,.317,.381,.413,.468,.457,.438,.384,.285,.295/ S=SSUNH(PHI,JD,RSTD) teta=90-wi(asin(s)) c write(*,*)teta QDOBSON=(1.-ALBEDO)*si0(jd)*cos(rad(teta))*(A(ICL)*S+B(ICL)*S*S) RETURN END REAL FUNCTION RAD(PHI) C ====================== C BOGENMASS IN WINKEL UMRECHNEN; REAL PHI RAD=PHI*3.141592654/180. RETURN END REAL FUNCTION WI(PHI) C ===================== REAL PHI WI=PHI*180./3.141592654 RETURN END SUBROUTINE SUNALT(T,RPHI,AL) C =========================== REAL ALRPHI INTEGER T C NOON SOLAR ALTITUDE (AL) AS FUNCTION OF DAY OF THE C YEAR (T) AND LATITUDE (RPHI), after SECKEL + BEAUDRY,1973 REAL A2,RT RT=T A1=RAD(23.87*SIN(RAD(RT-82.))) A2=RAD(RPHI) AL=(SIN(A2)*SIN(A1))+(COS(A2)*COS(A1)) AL=WI(ASIN(AL)) RETURN END C SUBROUTINE Q0SEBE(RPHI,T,QKW) C ============================= REAL QKW,RPHI INTEGER T c output: C daily means of SHORTWAVE RADIATION (DIRECT AND DIFFUSE) C WITHOUT CLOUDS C AFTER SECKEL + BEAUDRY, 1973 AND REED, 1977 c in W/m2 : QKW C as a function of c latitude : RPHI c Julian day : T c REAL RT,RADPHI,A0,A1,A2,B1,B2,Q RT=T Q=RAD(360.)/365.*(RT-21.) RADPHI=RAD(RPHI) IF (RPHI.LT.41.) THEN A0=-15.82+(326.87*COS(RADPHI)) A1=9.63+(192.44*COS(RAD(RPHI+90.))) A2=-0.64+(7.8*SIN(2.*RAD(RPHI-45.))) B1=-3.27+(108.7*SIN(RADPHI)) B2=-0.5+(14.42*COS(2.*RAD(RPHI-5.))) ELSE A0=342.61-(1.97*RPHI)-(0.018*RPHI*RPHI) A1=52.08-5.86*RPHI+(0.043*RPHI*RPHI) A2=-1.08-(0.47*RPHI)+(0.011*RPHI*RPHI) B1=-4.8+(2.46*RPHI)-(0.017*RPHI*RPHI) B2=-38.79+(2.43*RPHI)-(0.034*RPHI*RPHI) ENDIF QKW=(A0+A1*COS(Q)+A2*COS(2.*Q)+B1*SIN(Q)+B2*SIN(2.*Q)) C GLOBALSTRAHLUNG BEI WOLKENLOSEM HIMMEL RETURN END C FUNCTION QKCLOU(TC,MALT,MQKW) C ===================================== REAL TC,MALT,MQKW C CLOUD CORRECTION TO INSOLATION AFTER REED,1977 c input: c daily average of cloud cover in oktas : TC c noon solar altidue (degrees) : MALT c solar radiation under cloudless sky(W/m2): MQKW c c output: c solar radiation corrected due to clouds (W/m2) : QKCLOU C C EINGABE TC = TAGESMITTEL DER WOLKENBEDECKUNG IN ACHTEL C Z.B. 5/8 TC=5. C MALT= SONNENHOEHE ZU MITTAG IN GRAD C MQKW= GLOBALSTRAHLUNG OHNE WOLKEN c REAL TTC TTC=TC/(8.) IF (TTC.EQ.0.) THEN QKCLOU=MQKW RETURN ENDIF IF (TTC.LT..3) THEN QKCLOU=0.95*MQKW ELSE QKCLOU=(1.-(0.62*TTC)+(0.0019*MALT))*MQKW ENDIF RETURN END C FUNCTION JULDAY(MM,ID,IYYY) C *************************** PARAMETER (IGREG=15+31*(10+12*1582)) C AUS RECIPES, NA LOGO c c absolute Julian day c as a function of month (MM), day (ID) and year (IYYY) IF (IYYY.EQ.0) PAUSE 'There is no Year Zero.' IF (IYYY.LT.0) IYYY=IYYY+1 IF (MM.GT.2) THEN JY=IYYY JM=MM+1 ELSE JY=IYYY-1 JM=MM+13 ENDIF JULDAY=INT(365.25*JY)+INT(30.6001*JM)+ID+1720995 IF (ID+31*(MM+12*IYYY).GE.IGREG) THEN JA=INT(0.01*JY) JULDAY=JULDAY+2-JA+INT(0.25*JA) ENDIF RETURN END SUBROUTINE JULD90(M,ID,ITAG90) C ****************************** c c Julian day of the year (1-365) (ITAG90) as a function of c month (M) and day of the month (ID) C JULIANISCHER TAG IM JAHR 1990 und in jedem anderen Jahr C VON MIR JAHR=1990 ITAG90=JULDAY(M,ID,JAHR)-JULDAY(12,31,JAHR-1) RETURN END SUBROUTINE HAABB(XPHI,JULTAG,AA,BB,ET) C =================================== INTEGER JULTAG REAL AA,BB,XPHI,ET C VALUES AA AND BB IN FORMULAR FOR SIN ALPHA AS FUNCTION C OF LATITUDE AND JULIAN DAY, AFTER DAVIES FROM SMITH & DOBSON C REAL RPHI,THETA,PI,DEC,TH2,TH3,CT,ST,C2T,S2T,C3T,S3T RPHI=RAD(XPHI) PI=3.141593 THETA=2.*PI*REAL(JULTAG-1)/365. TH2=2.*THETA TH3=3.*THETA CT=COS(THETA) ST=SIN(THETA) C2T=COS(TH2) S2T=SIN(TH2) C3T=COS(TH3) S3T=SIN(TH3) DEC=0.006918-0.399912*CT+0.070257*ST-0.006759*C2T+0.000907*S2T- *0.002697*C3T+0.001480*S3T AA=SIN(DEC)*SIN(RPHI) BB=COS(DEC)*COS(RPHI) ET=0.000075+0.001868*CT-0.032077*ST-0.014615*C2T-0.040840*S2T RETURN END C FUNCTION SSUNH(PHI,JD,RSTD) C *************************** INTEGER JD REAL PHI,RSTD C c SSUNH is the sine of the solar elevation angle as a function of c latitude (PHI) c Julian day (JD) c hour of the day (RSTD) C SSUNH IST DER SINUS DES SONNENHOEHENWINKELS AM C JULIANISCHEN TAG JD ZUR STUNDE ISTD IN DER GEOGR. BREITE PHI C REAL RST,AA,BB,ET,HA HA=ABS(12.-RSTD)*3.14159/12. CALL HAABB(PHI,JD,AA,BB,ET) SSUNH=AA+BB*COS(HA) RETURN END SUBROUTINE SAUFUN(JD,PHI,SAUF,SUNT) C =================================== REAL PHI,SAUF,SUNT INTEGER JD C c input: c Julian day (JD) c latitude (PHI) c c output: c local time of sunrise (SAUF) and sunset (SUNT) c C BERECHNET FUER DEN JULIANISCHEN TAG JD IN DER GEOGRAPHISCHEN C BREITE PHI DIE UHRZEITEN DES SONNENAUF- UND -UNTERGANGS C C Z.ZT. WERDEN SCHWEINEREIEN POLWAERTS DER C POLARKREISE NICHT BERUECKSICHTIGT C C LAENGENKORREKTUREN INNERHALB VON 1 STUNDEN ZEITZONEN C MUESSTEN ADDIERT WERDEN: C HIER STEHT DIE SONNE MITTAGS (12 UHR) AM HOECHSTEN PI=3.14159 CALL HAABB(PHI,JD,AA,BB,ET) CC=-AA/BB H1=ACOS(CC) HH=H1*12./PI HET=12.-ET SAUF=ABS(HET-HH) SUNT=ABS(HET+HH) 2 CONTINUE RETURN END C SUBROUTINE CLOUDA(JD,PHI,NCL,XM1CL,XM2CL,I1,I2,IG) C *********************************************** INTEGER JD,I1,I2 REAL PHI,XM1CL,XM2CL INTEGER NCL(0:23) C C DIE FUER DIE PARAMETRISIERUNG DER KURZWELLIGEN EINSTRAHLUNG C NACH REED RELEVANTE MITTLERE WOLKENBEDECKUNG EINES TAGES C WIRD AUS DEN STUENDLICHEN C BEOBACHTUNGEN DER GESAMTBEDECKUNG ERMITTELT. C C EINGABE: C JD : JULIANISCHER TAG IM JAHR C PHI : GEOGR.BREITE C NCL : ENTHAELT DIE STUENDLICHEN WOLKENBEOBACHTUNGEN C ALS INTEGERZAHLEN, 0 BIS 8, 9 = VERDECKT, 99 = NICHT C VORHANDEN ODER AEHNLICH C C AUSGABE: C C I1,I2: ERSTE UND LETZTE VOLLE STUNDE MIT SONNE UEBER DEM HORIZONT C D.H. SIN ALPHA > 0.1 C XM1CL, XM2CL: MITTELWERT UND STANDARDABWEICHUNG AUS DEN C BEOBACHTUNGEN VON I1 BIS I2 C C C IG : GESAMTANZAHL AUSGEWERTETER BEOBACHTUNGEN REAL X(24) C CALL SAUFUN(JD,PHI,SAUF,SUNT) I1=0 I2=0 C FINDE ERSTE VOLLE STUNDE MIT SIN ALPHA > 0.1 DO 8 IH=1,12 RH=IH SA=SSUNH(PHI,JD,RH) IF((I1.EQ.0).AND.(SA.GT.0.1))THEN I1=IH GOTO 9 ENDIF 8 CONTINUE 9 CONTINUE C FINDE LETZTE VOLLE STUNDE MIT SIN ALPHA > 0.1 DO 18 IH=13,24 RH=IH SA=SSUNH(PHI,JD,RH) IF((I2.EQ.0).AND.(SA.LT.0.1))THEN I2=IH-1 GOTO 19 ENDIF 18 CONTINUE 19 CONTINUE C I1=NINT(ANINT(SAUF+0.5)) C I2=NINT(ANINT(SUNT-0.5)) IZ=0 DO 10 I=I1,I2 IF((NCL(I).GE.0).AND.(NCL(I).LE.8))THEN IZ=IZ+1 X(IZ)=REAL(NCL(I)) ENDIF 10 CONTINUE IF (IZ.EQ.0)THEN XM1CL=0. XM2CL=0. IG=0 RETURN ENDIF IF (IZ.LT.2)THEN XM1CL=X(1) XM2CL=0. IG=1 RETURN ENDIF CALL MOMENT(X,IZ,XM1CL,A1,XM2CL,A2,A3,A4) IG=IZ RETURN END C SUBROUTINE QDODAY(ALB,JD,PHI,NCL,XMCT,XSCT,QDOB,XMQ,XSQ,I1,I2,IG) C ************************************************************* INTEGER JD,I1,I2 REAL PHI,XMCT,XSCT,ALB,QDOB(0:23),XMQ,XSQ INTEGER NCL(0:23) C C DIE FUER DIE PARAMETRISIERUNG DER KURZWELLIGEN EINSTRAHLUNG C NACH REED RELEVANTE MITTLERE WOLKENBEDECKUNG EINES TAGES C WIRD AUS DEN STUENDLICHEN C BEOBACHTUNGEN DER GESAMTBEDECKUNG ERMITTELT. C C FUER JEDE STUNDE MIT SIN ALPHA > 0.1 WIRD DIE ABSORBIERTE C STRAHLUNG NACH DEM OKTA MODELL VON DOBSON UND SMITH BERECHNET C UND IN DEM ARRAY QDOB UEBERGEBEN. MITTELWERT UND STANDARDABWEICHUNG C FUER DEN TAG WERDEN IN XMQ UND XSQ UEBERGEBEN. C EINGABE: C JD : JULIANISCHER TAG IM JAHR C PHI : GEOGR.BREITE C NCL : ENTHAELT DIE STUENDLICHEN WOLKENBEOBACHTUNGEN C ALS INTEGERZAHLEN, 0 BIS 8, 9 = VERDECKT, 99 = NICHT C VORHANDEN ODER AEHNLICH C C AUSGABE: C C I1,I2: ERSTE UND LETZTE VOLLE STUNDE MIT SONNE UEBER DEM HORIZONT C D.H. SIN ALPHA > 0.1 C XMCT, XSCT: MITTELWERT UND STANDARDABWEICHUNG AUS DEN C WOLKENBEOBACHTUNGEN VON I1 BIS I2 C C IG : GESAMTANZAHL AUSGEWERTETER BEOBACHTUNGEN C C QDOB: ARRAY MIT DEN VOM OZEAN ABSORBIERTEN STRAHLUNGSFLUESSEN C FUER DIE STUNDEN, WO SIN ALPHA > 0.1 IST, SONST C STEHT EINE NULL ODER -99 WENN KEINE WOLKEN DA. C XMQ,XSQ : MITTELWERT UND STANDARDABWEICHUNG DER STUNDENWERTE C DER STRAHLUNG FUER I1 BIS I2. REAL X(24),Y(24) DO 99 I=0,23 QDOB(I)=0. 99 CONTINUE C CALL SAUFUN(JD,PHI,SAUF,SUNT) I1=0 I2=0 C FINDE ERSTE VOLLE STUNDE MIT SIN ALPHA > 0.1 DO 8 IH=1,12 RH=IH SA=SSUNH(PHI,JD,RH) IF((I1.EQ.0).AND.(SA.GT.0.1))THEN I1=IH GOTO 9 ENDIF 8 CONTINUE 9 CONTINUE C FINDE LETZTE VOLLE STUNDE MIT SIN ALPHA > 0.1 DO 18 IH=13,24 RH=IH SA=SSUNH(PHI,JD,RH) IF((I2.EQ.0).AND.(SA.LT.0.1))THEN I2=IH-1 GOTO 19 ENDIF 18 CONTINUE 19 CONTINUE C I1=NINT(ANINT(SAUF+0.5)) C I2=NINT(ANINT(SUNT-0.5)) IZ=0 DO 10 I=I1,I2 RI=I IF((NCL(I).GE.0).AND.(NCL(I).LE.8))THEN IZ=IZ+1 X(IZ)=REAL(NCL(I)) QDOB(I)=QDOBSON(PHI,JD,RI,NCL(I),ALB) Y(IZ)=QDOB(I) ELSE QDOB(I)=-99. ENDIF 10 CONTINUE IF (IZ.EQ.0)THEN XMCT=0. XSCT=0. IG=0 XMQ=0. XSQ=0. RETURN ENDIF IF (IZ.LT.2)THEN XMCT=X(1) XSCT=0. XMQ=Y(1) XSQ=0. IG=1 RETURN ENDIF CALL MOMENT(X,IZ,XMCT,A1,XSCT,A2,A3,A4) CALL MOMENT(Y,IZ,XMQ,A1,XSQ,A2,A3,A4) IG=IZ RETURN END C REAL FUNCTION QTAG(PHI,JD) C************************************************************** C QTAG IS THE DAILY SOLAR RADIATION AT THE TOP OF THE C ATMOSPHERE AS A FUNCTION OF C LATITUDE : PHI C JULIAN DAY : JD C EXPRESSED IN WM**-2 AVERAGED OVER 24 HOURS C SEE E.G. THE TEXTBOOK BY KONDRATYEV, P.309 C************************************************************** REAL PHI INTEGER JD BPHI=RAD(PHI) D=DEKLIJ(JD) BDEK=RAD(D) H0=HAUF(PHI,D) BH0=RAD(H0) S=SI0(JD) ERG1=(BH0*SIN(BPHI)*SIN(BDEK)) + (COS(BPHI)*COS(BDEK)*SIN(BH0)) QTAG=ERG1/3.1415927*SI0(JD) END C REAL FUNCTION HAUF(PHI,DEKLI) C****************************************************************** C CALCULATES THE HOUR ANGLE OF SUNRISE AND SUNSET (THETA=90) C IN DEGREES C SEE FOR EXAMPLE THE TEXTBOOK OF KONDRATYEV, P.308 C INPUT C LATITUDE : PHI C DEKLINATION OF THE SUN : DEKLI C******************************************************************** REAL PHI,DEKLI BPHI=RAD(PHI) BDEK=RAD(DEKLI) ERG1=TAN(BPHI)*TAN(BDEK) HAUF=WI(ACOS(-ERG1)) END SUBROUTINE MOMENT(DATA,N,AVE,ADEV,SDEV,VAR,SKEW,CURT) DIMENSION DATA(N) IF(N.LE.1)PAUSE 'N must be at least 2' S=0. DO 11 J=1,N S=S+DATA(J) 11 CONTINUE AVE=S/N ADEV=0. VAR=0. SKEW=0. CURT=0. DO 12 J=1,N S=DATA(J)-AVE ADEV=ADEV+ABS(S) P=S*S VAR=VAR+P P=P*S SKEW=SKEW+P P=P*S CURT=CURT+P 12 CONTINUE ADEV=ADEV/N VAR=VAR/(N-1) SDEV=SQRT(VAR) IF(VAR.NE.0.)THEN SKEW=SKEW/(N*SDEV**3) CURT=CURT/(N*VAR**2)-3. ELSE PAUSE 'no skew or kurtosis when zero variance' ENDIF RETURN END SUBROUTINE INTERP (X,Y,F2,F1,M,N,A,RES) DIMENSION A(M,N),F1(N),F2(M) K1=M-1 IF(Y.EQ.F2(M)) M1=M DO 3 I=1,K1 IF(Y.GE.F2(I).AND.Y.LT.F2(I+1)) GOTO 4 3 CONTINUE 4 M1=I K2=N-1 IF(X.EQ.F1(N)) N1=N DO 5 I=1,K2 IF(X.GE.F1(I).AND.X.LT.F1(I+1)) GOTO 6 5 CONTINUE 6 N1=I R1=A(M1,N1)+(A(M1,N1+1)-A(M1,N1))/ /(F1(N1+1)-F1(N1))*(X-F1(N1)) R2=A(M1+1,N1)+(A(M1+1,N1+1)-A(M1+1,N1))/ /(F1(N1+1)-F1(N1))*(X-F1(N1)) RES=R1+(R2-R1)/(F2(M1+1)-F2(M1))*(Y-F2(M1)) RETURN END C******************************************************************** C s h o r t w a v e r a d i a t i o n C monthly mean values C Sergey Malevsky scheme 1991 -version C Subroutine prepared and tested in January 12-18,1992 St. Petersburg C and in March 13-19,1992 Moscow C by Sergey Malevsky and Sergey Gulev C******************************************************************** C AM - month REAL " C SR - latitude REAL " C TA - air temperature deg,C REAL " C EZ - air humidity mb REAL " C ANN- low cloudness percent/10 REAL " C ANO- common cloudness percent/10 REAL " C QS - total radiation W/M**2, REAL output C RPOGL-accumulate radiation W/M**2, REAL " C IE - =0, if two levels cloudness INTEGER input C =1, if only common C******************************************************************* SUBROUTINE RSWM (AM,SR,TA,EZ,ANN,ANO,QS,RPOGL,IE) DIMENSION HSOL(10),ANR1(14),ALB(10,14),FI(71),QQ0(12,71), *TAR(12),FN(14,14),FN1(14,14,10),FN0(14,10),q1(71), *qq01(12,36),qq02(12,35), *alfa(10,26),pal(26) DATA PAL/0.60,0.61,0.62,0.63,0.64,0.65,0.66,0.67,0.68,0.69, * 0.70,0.71,0.72,0.73,0.74,0.75,0.76,0.77,0.78,0.79, * 0.80,0.81,0.82,0.83,0.84,0.85/ DATA HSOL/5.,10.,20.,30.,40.,50.,60.,70.,80.,90./ DATA ANR1/0.,1.,2.,3.,4.,5.,6.,7.,7.5,8.,8.5,9.,9.5,10./ DATA TAR/-20.,-15.,-10.,-5.,0.,5.,10.,15.,20.,25.,30.,35./ DATA ALB/0.19,0.21,0.16,0.11,0.08,0.07,0.06,0.06,0.05,0.05, * 0.19,0.21,0.16,0.11,0.08,0.07,0.06,0.06,0.05,0.05, * 0.19,0.21,0.16,0.11,0.08,0.07,0.06,0.06,0.05,0.05, * 0.18,0.20,0.16,0.11,0.08,0.07,0.06,0.06,0.05,0.05, * 0.17,0.19,0.16,0.11,0.08,0.07,0.06,0.06,0.05,0.05, * 0.16,0.18,0.16,0.11,0.08,0.07,0.06,0.06,0.05,0.05, * 0.16,0.17,0.15,0.11,0.09,0.08,0.07,0.07,0.06,0.06, * 0.14,0.15,0.14,0.11,0.09,0.08,0.07,0.07,0.06,0.06, * 0.13,0.14,0.13,0.11,0.09,0.08,0.07,0.07,0.06,0.06, * 0.12,0.13,0.12,0.11,0.09,0.08,0.07,0.07,0.06,0.06, * 0.12,0.12,0.11,0.11,0.09,0.08,0.07,0.07,0.06,0.06, * 0.11,0.11,0.11,0.10,0.09,0.08,0.07,0.07,0.06,0.06, * 0.10,0.10,0.10,0.10,0.09,0.08,0.07,0.07,0.06,0.06, * 0.10,0.10,0.10,0.10,0.09,0.08,0.07,0.07,0.06,0.06/ DATA FN0/1.00,1.02,1.01,0.97,0.93,0.84,0.76, * 0.67,0.63,0.59,0.56,0.52,0.46,0.39, * 1.00,1.02,1.01,0.97,0.93,0.84,0.76, * 0.67,0.63,0.59,0.55,0.51,0.44,0.37, * 1.00,1.02,1.01,0.97,0.92,0.84,0.76, * 0.67,0.63,0.58,0.54,0.49,0.43,0.35, * 1.00,1.01,1.00,0.97,0.92,0.84,0.77, * 0.68,0.64,0.59,0.55,0.50,0.44,0.37, * 1.00,1.00,0.99,0.96,0.92,0.85,0.78, * 0.70,0.65,0.61,0.56,0.51,0.46,0.40, * 1.00,1.00,0.98,0.96,0.92,0.85,0.79, * 0.71,0.67,0.63,0.58,0.53,0.48,0.42, * 1.00,0.99,0.97,0.95,0.92,0.85,0.79, * 0.73,0.68,0.64,0.60,0.55,0.49,0.44, * 1.00,0.99,0.97,0.95,0.92,0.85,0.80, * 0.74,0.70,0.66,0.61,0.57,0.51,0.46, * 1.00,0.98,0.96,0.94,0.91,0.86,0.80, * 0.75,0.71,0.67,0.63,0.59,0.53,0.48, * 1.00,0.98,0.96,0.94,0.91,0.87,0.81, * 0.76,0.72,0.68,0.65,0.61,0.56,0.50/ DATA ALFA/0.64,0.68,0.74,0.78,0.81,0.84,0.85,0.85,0.86,0.86, * 0.66,0.70,0.76,0.79,0.82,0.85,0.86,0.86,0.87,0.87, * 0.68,0.72,0.77,0.81,0.84,0.86,0.87,0.87,0.87,0.87, * 0.70,0.74,0.79,0.82,0.85,0.87,0.88,0.88,0.88,0.88, * 0.72,0.76,0.80,0.83,0.86,0.88,0.89,0.89,0.89,0.89, * 0.74,0.78,0.82,0.85,0.87,0.89,0.90,0.90,0.91,0.91, * 0.76,0.80,0.84,0.86,0.88,0.90,0.91,0.91,0.92,0.92, * 0.78,0.82,0.85,0.88,0.90,0.91,0.92,0.92,0.93,0.93, * 0.80,0.84,0.87,0.89,0.91,0.92,0.93,0.93,0.93,0.93, * 0.82,0.86,0.88,0.90,0.92,0.93,0.94,0.94,0.94,0.94, * 0.84,0.88,0.90,0.92,0.93,0.94,0.95,0.95,0.95,0.95, * 0.88,0.90,0.92,0.94,0.94,0.95,0.96,0.96,0.96,0.96, * 0.91,0.93,0.94,0.95,0.96,0.96,0.97,0.97,0.97,0.97, * 0.94,0.95,0.96,0.97,0.97,0.98,0.98,0.98,0.98,0.98, * 0.97,0.98,0.98,0.98,0.99,0.99,0.99,0.99,0.99,0.99, * 1.00,1.00,1.00,1.00,1.00,1.00,1.00,1.00,1.00,1.00, * 1.03,1.03,1.02,1.02,1.02,1.02,1.02,1.01,1.01,1.01, * 1.06,1.05,1.04,1.04,1.03,1.03,1.03,1.03,1.03,1.03, * 1.09,1.08,1.07,1.06,1.05,1.05,1.05,1.05,1.05,1.05, * 1.12,1.10,1.09,1.08,1.07,1.06,1.06,1.06,1.06,1.06, * 1.15,1.13,1.11,1.10,1.09,1.08,1.08,1.08,1.08,1.08, * 1.16,1.14,1.13,1.12,1.11,1.10,1.10,1.10,1.10,1.10, * 1.17,1.16,1.15,1.14,1.13,1.12,1.12,1.12,1.12,1.12, * 1.18,1.18,1.17,1.17,1.16,1.15,1.14,1.13,1.13,1.13, * 1.19,1.19,1.19,1.19,1.18,1.17,1.16,1.15,1.15,1.15, * 1.20,1.21,1.21,1.21,1.20,1.19,1.18,1.17,1.17,1.17/ data QQ01/ *29.5,22.4,14.3, 6.9, 2.7, 1.2, 1.7, 4.6,10.6,18.6,26.8,31.4, *29.8,23.0,15.1, 7.8, 3.4, 1.7, 2.3, 5.4,11.5,19.4,27.2,31.5, *30.1,23.6,16.0, 8.6, 4.1, 2.4, 3.0, 6.2,12.4,20.1,27.6,31.7, *30.4,24.2,16.8, 9.5, 4.9, 3.0, 3.7, 7.0,13.2,20.8,28.0,31.9, *30.6,24.8,17.6,10.4, 5.7, 3.7, 4.4, 7.8,14.0,21.4,28.4,32.1, *30.9,25.3,18.3,11.2, 6.5, 4.5, 5.2, 8.7,14.9,22.1,28.8,32.2, *31.6,25.8,19.1,12.1, 7.3, 5.2, 6.0, 9.5,15.7,22.7,29.1,32.4, *31.3,26.3,19.8,13.0, 8.2, 6.0, 6.8,10.4,16.5,23.3,29.4,32.9, *31.5,26.8,20.5,13.8, 9.0, 6.8, 7.6,11.2,17.3,23.9,29.7,32.6, *31.6,27.2,21.2,14.6, 9.9, 7.7, 8.4,12.1,18.0,24.4,29.9,32.6, *31.7,27.6,21.9,15.5,10.7, 8.5, 9.3,12.9,18.8,24.9,30.1,32.7, *31.8,27.9,22.5,16.3,11.6, 9.4,10.1,13.8,29.5,25.4,30.3,32.6, *31.9,28.2,23.1,17.1,12.4,10.2,11.0,14.6,20.2,25.8,30.4,32.6, *31.9,28.5,23.6,17.8,13.3,11.1,11.9,15.4,20.8,26.2,30.5,32.5, *31.8,28.8,24.2,18.6,14.1,12.0,12.7,16.2,21.5,26.6,30.6,32.4, *31.8,29.0,24.7,19.3,15.0,12.8,13.6,17.0,22.1,27.0,30.6,32.3, *31.7,29.2,25.2,20.1,15.8,13.7,14.4,17.8,22.7,27.3,30.6,32.1, *31.5,29.3,25.6,20.8,16.6,14.6,15.3,18.6,23.2,27.5,30.6,31.9, *31.4,29.4,26.0,21.4,17.4,15.4,16.1,19.3,23.8,27.8,30.5,31.7, *31.1,29.5,26.4,22.1,18.2,16.2,16.9,20.0,24.2,28.0,30.4,31.4, *30.9,29.5,26.7,22.7,19.0,17.1,17.7,20.7,24.7,28.1,30.2,31.1, *30.6,29.5,27.0,23.3,19.8,17.9,18.5,21.4,25.1,28.2,30.1,30.7, *30.3,29.4,27.3,23.8,20.5,18.7,19.3,22.0,25.5,28.3,29.8,30.3, *29.9,29.4,27.5,24.3,21.2,19.5,20.0,22.6,25.9,28.4,29.6,29.9, *29.5,29.8,27.7,24.8,21.9,20.2,20.8,23.2,26.2,28.4,29.3,29.4, *29.1,29.1,27.8,25.3,22.6,21.0,21.5,23.8,26.5,28.3,28.9,28.9, *28.6,28.8,27.9,25.7,23.2,21.7,22.2,24.3,26.8,28.3,28.6,28.4, *28.2,28.6,28.0,26.1,23.8,22.4,22.9,24.8,27.0,28.1,28.0,27.9, *27.6,28.3,28.0,26.5,24.4,23.1,23.5,25.3,27.2,28.0,27.7,27.3, *27.1,28.0,28.0,26.8,25.0,23.8,24.1,25.7,27.3,27.8,27.2,26.6, *26.5,27.7,28.0,27.1,25.5,24.4,24.7,26.1,27.4,27.6,26.7,26.0, *25.8,27.3,27.9,27.4,26.0,25.0,25.3,26.5,27.5,27.3,26.2,25.3, *25.2,26.8,27.8,27.6,26.5,25.6,25.8,26.8,27.5,27.0,25.6,24.6, *24.5,26.4,27.6,27.8,26.9,26.2,26.3,27.2,27.5,26.7,25.0,23.8, *23.8,25.9,27.4,27.9,27.3,26.7,26.8,27.4,27.5,26.3,24.4,23.1, *23.1,25.4,27.2,28.0,27.7,27.2,27.3,27.7,27.4,25.9,23.7,22.4/ data QQ02/ *22.3,24.8,26.9,28.1,28.0,27.6,27.7,27.9,27.3,25.5,23.0,21.6, *21.5,24.2,26.6,28.1,28.3,28.1,28.0,28.0,27.1,25.0,22.3,20.7, *20.7,23.6,26.3,28.1,28.6,28.5,28.4,28.2,26.9,24.5,21.6,19.9, *19.9,23.0,25.9,28.1,28.8,28.8,28.7,28.2,26.7,24.0,20.8,19.0, *19.1,22.3,25.5,28.0,29.0,29.2,29.0,28.3,26.4,23.4,20.0,18.2, *18.2,21.6,25.0,27.9,29.2,29.5,29.2,28.3,26.1,22.8,19.2,17.3, *17.4,20.8,24.5,27.7,29.3,29.7,29.4,28.3,25.8,22.2,18.4,16.4, *16.5,20.1,24.0,27.6,29.4,30.0,29.6,28.2,25.4,21.6,17.6,15.5, *15.6,19.3,23.5,27.3,29.4,30.2,29.8,28.1,25.0,20.9,16.7,14.6, *14.7,18.5,22.9,27.1,29.4,30.3,29.9,28.0,24.6,20.2,15.8,13.7, *13.8,17.7,22.3,26.8,29.4,30.5,30.0,27.9,24.1,19.4,15.0,12.7, *12.9,16.9,21.7,26.4,29.4,30.6,30.0,27.7,23.6,18.7,14.1,11.8, *11.9,16.1,21.0,26.1,29.3,30.6,30.0,27.4,23.0,17.9,13.2,10.9, *11.0,15.2,20.3,25.7,29.2,30.7,30.0,27.2,22.5,17.1,12.3,10.0, *10.1,14.3,19.6,25.2,29.0,30.7,29.9,26.9,21.9,16.3,11.4, 9.1, * 9.2,13.4,18.9,24.8,28.8,30.6,29.8,26.5,21.3,15.5,10.5, 8.2, * 8.3,12.6,18.1,24.3,28.6,30.6,29.7,26.2,20.6,14.7, 9.6, 7.3, * 7.4,11.7,17.3,23.7,28.4,30.5,29.6,25.8,19.9,13.8, 8.7, 6.4, * 6.6,10.8,16.5,23.2,28.1,30.4,29.4,25.4,19.2,13.0, 7.8, 5.6, * 5.7, 9.9,15.7,22.6,27.8,30.3,29.2,24.9,18.5,12.1, 7.0, 4.8, * 4.9, 9.0,14.9,22.0,27.4,30.1,29.0,24.4,17.8,11.2, 6.1, 4.0, * 4.1, 8.1,14.0,21.4,27.1,30.0,28.8,23.9,17.0,10.3, 5.1, 3.2, * 3.3, 7.2,13.2,20.7,26.8,29.8,28.5,23.4,16.2, 9.4, 4.5, 2.5, * 2.6, 6.3,12.3,20.0,26.4,29.6,28.2,22.9,15.4, 8.6, 3.7, 1.8, * 1.9, 5.5,11.4,19.3,26.0,29.4,28.0,22.3,14.6, 7.7, 2.9, 1.3, * 1.3, 4.6,10.5,18.6,25.6,29.3,27.7,21.8,13.8, 6.8, 2.2, .8, * 0.8, 3.8, 9.6,17.9,25.9,29.2,27.5,21.2,12.9, 6.0, 1.6, .4, * 0.4, 3.1, 8.7,17.2,24.8,29.2,27.4,20.6,12.1, 5.1, 1.0, .1, * 0.1, 2.3, 7.8,16.4,24.5,29.3,27.3,20.0,11.2, 4.3, .6, .0, * 0.0, 1.7, 6.9,15.7,24.3,29.5,27.4,19.5,10.4, 3.5, .2, .0, * 0.0, 1.1, 6.0,15.0,24.2,29.7,27.5,19.0, 8.5, 2.7, .0, .0, * 0.0, 0.6, 5.2,14.2,24.3,30.0,27.6,18.5, 8.6, 2.0, .0, .0, * 0.0, 0.2, 4.3,13.6,24.3,30.1,27.8,18.3, 7.8, 1.3, .0, .0, * 0.0, 0.0, 3.4,13.0,24.4,30.3,27.9,18.2, 6.9, .8, .0, .0, * 0.0, 0.0, 2.6,12.5,24.4,30.4,28.0,18.1, 6.1, .3, .0, .0/ do 108 j=1,36 do 108 i=1,12 108 qq0(i,j)=qq01(i,j) do 109 j=37,71 j2=j-36 do 109 i=1,12 109 qq0(i,j)=qq02(i,j2) call WIS (am,sr,H) IF(H.LE.5.) H=5.01 c********************************************************************* c P2 - constant calculations c 1) P2=f(ez) 2) P2=f(ta) c choose your situation and comment or uncomment others c********************************************************************* c 1) P2=f(ez) - for North Atlantic c********************************************************************* ccc p2=0.829-0.0078*ez+0.000115*ez*ez c********************************************************************* c 1) P2=f(ez) - for South Atlantic, Indian and Pacific c********************************************************************* ccc p2=0.797-0.0032*ez+0.000034*ez*ez c********************************************************************* c 2) P2=f(ta) - for North Atlantic c********************************************************************* p2=0.799-0.0037*ta c********************************************************************* c 2) P2=f(ta) - for South Atlantic, Indian and Pacific c********************************************************************* ccc p2=0.785-0.0018*ta c********************************************************************* c 2) P2=f(ta) - for World Ocean c********************************************************************* ccc p2=0.790-0.0030*ta c********************************************************************* c possible clear sky incoming radiation under P2=0.75 c********************************************************************* do 199 if=1,71 199 fi(if)=-60.+2.*float(if-1) im=int(am) do 299 if=1,71 299 q1(if)=qq0(im,if) do 399 if=1,70 if (fi (if).lt.sr.and.fi (if+1).ge.sr) go to 499 go to 399 499 Q00=q1(if)+((sr-fi (if))/(fi (if+1)-fi (if)))*(q1(if+1)-q1(if)) 399 continue c********************************************************************* c clear sky incoming radiation under real P2 c********************************************************************* H1=H CALL INTERP(P2,H1,HSOL,PAL ,10,26,ALFA,AAA) Q0=AAA*Q00 c********************************************************************* ifl=0 IF(H1.LE.5.01) ifl=1 c********************************************************************* c albedo calculation -- albedo = F (solar high, total cloud cover) c********************************************************************* CALL INTERP(ANO,H1,HSOL,ANR1,10,14,ALB,ALP) c********************************************************************* c choice the way (TOTAL/TOTAL+LOW) according IE-parameter c********************************************************************* IF(IE.EQ.1) GOTO 16 c********************************************************************* c 1-st way --> T O T A L + L O W c********************************************************************* c preparing (overturning) of 2-levels cloudness matrix (FN1) c********************************************************************* CALL DOFM (FN1) c********************************************************************* c cloud cover factor (RF) determination by interpolation c********************************************************************* DO 10 I=1,9 IF(H1.GT.HSOL(I).AND.H1.LE.HSOL(I+1)) GOTO 11 10 CONTINUE 11 N=I DO 12 L=1,14 DO 12 J=1,14 12 FN(L,J)=FN1(J,L,N) CALL INTERP (ANN,ANO,ANR1,ANR1,14,14,FN,S1) N=N+1 DO 13 L=1,14 DO 13 J=1,14 13 FN(L,J)=FN1(J,L,N) CALL INTERP (ANN,ANO,ANR1,ANR1,14,14,FN,S2) RF=S1+(S2-S1)/(HSOL(N)-HSOL(N-1))*(H1-HSOL(N-1)) GOTO 17 c********************************************************************* c 2-nd way --> T O T A L o n l y c********************************************************************* 16 CONTINUE c********************************************************************* c cloud cover factor (RF) determination by interpolation c********************************************************************* CALL INTERP(H1,ANO,ANR1,HSOL,14,10,FN0,RF) 17 CONTINUE c********************************************************************* c incoming shortwave radiation c********************************************************************* QS=Q0*RF*100./8.64 c********************************************************************* c absorbed shortwave radiation c********************************************************************* RPOGL=QS*(1-ALP) c********************************************************************* c small solar high c********************************************************************* if(ifl.eq.0) go to 6789 qs=0. rbl=0. rpogl=0. 6789 continue RETURN END c********************************************************************* c forming of two-level cloud cover coefficints c preparing (overturning) of 2-levels cloudness matrix (F), c including problem of zero-values influence on interpolation c********************************************************************* SUBROUTINE DOFM(F) DIMENSION F(14,14,10),D1(14,14),D2(14,14),D3(14,14),D4(14,14), *D5(14,14),D6(14,14),D7(14,14),D8(14,14),D9(14,14),D10(14,14) DATA D1 */1.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 1.02,1.02,0.00,0.00,0.00,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.99,1.00,1.01,0.00,0.00,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.97,0.97,0.96,0.96,0.00,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.93,0.93,0.92,0.91,0.89,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.87,0.87,0.86,0.84,0.81,0.79,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.82,0.82,0.81,0.78,0.76,0.73,0.70, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.76,0.76,0.75,0.72,0.70,0.68,0.65, * 0.63,0.00,0.00,0.00,0.00,0.00,0.00, * 0.73,0.73,0.71,0.69,0.66,0.64,0.62, * 0.59,0.58,0.00,0.00,0.00,0.00,0.00, * 0.70,0.70,0.68,0.66,0.63,0.60,0.58, * 0.55,0.54,0.53,0.00,0.00,0.00,0.00, * 0.66,0.66,0.64,0.62,0.59,0.57,0.55, * 0.52,0.51,0.50,0.49,0.00,0.00,0.00, * 0.61,0.61,0.60,0.58,0.55,0.53,0.51, * 0.48,0.47,0.46,0.45,0.44,0.00,0.00, * 0.58,0.58,0.56,0.54,0.52,0.49,0.47, * 0.45,0.44,0.43,0.42,0.41,0.40,0.00, * 0.54,0.54,0.52,0.50,0.48,0.45,0.43, * 0.41,0.40,0.39,0.38,0.37,0.36,0.35/ DATA D2 */1.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 1.02,1.02,0.00,0.00,0.00,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.99,1.00,1.01,0.00,0.00,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.97,0.97,0.96,0.96,0.00,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.93,0.93,0.92,0.92,0.90,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.88,0.88,0.87,0.85,0.82,0.80,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.82,0.82,0.81,0.78,0.76,0.73,0.70, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.76,0.76,0.75,0.72,0.70,0.68,0.65, * 0.63,0.00,0.00,0.00,0.00,0.00,0.00, * 0.73,0.73,0.72,0.69,0.66,0.64,0.62, * 0.59,0.58,0.00,0.00,0.00,0.00,0.00, * 0.70,0.70,0.69,0.66,0.63,0.60,0.58, * 0.55,0.54,0.53,0.00,0.00,0.00,0.00, * 0.65,0.65,0.64,0.61,0.58,0.56,0.54, * 0.51,0.50,0.49,0.48,0.00,0.00,0.00, * 0.60,0.60,0.59,0.57,0.54,0.52,0.50, * 0.47,0.46,0.45,0.44,0.43,0.00,0.00, * 0.57,0.57,0.55,0.53,0.51,0.48,0.46, * 0.44,0.42,0.41,0.40,0.39,0.38,0.00, * 0.53,0.53,0.51,0.49,0.47,0.44,0.42, * 0.40,0.39,0.38,0.37,0.36,0.35,0.34/ DATA D3 */1.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 1.01,1.02,0.00,0.00,0.00,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.99,1.00,1.01,0.00,0.00,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.97,0.97,0.96,0.96,0.00,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.93,0.93,0.92,0.91,0.90,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.89,0.89,0.88,0.86,0.83,0.80,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.83,0.83,0.82,0.79,0.76,0.74,0.71, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.77,0.77,0.76,0.73,0.71,0.69,0.66, * 0.64,0.00,0.00,0.00,0.00,0.00,0.00, * 0.73,0.73,0.72,0.70,0.67,0.64,0.62, * 0.59,0.58,0.00,0.00,0.00,0.00,0.00, * 0.70,0.70,0.68,0.66,0.63,0.60,0.57, * 0.54,0.53,0.51,0.00,0.00,0.00,0.00, * 0.65,0.65,0.64,0.62,0.58,0.56,0.53, * 0.50,0.49,0.48,0.47,0.00,0.00,0.00, * 0.60,0.60,0.59,0.56,0.54,0.51,0.49, * 0.46,0.45,0.44,0.42,0.41,0.00,0.00, * 0.56,0.56,0.54,0.52,0.50,0.48,0.46, * 0.43,0.42,0.40,0.39,0.38,0.37,0.00, * 0.53,0.53,0.50,0.48,0.46,0.44,0.42, * 0.40,0.39,0.37,0.36,0.35,0.34,0.33/ DATA D4 */1.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 1.01,1.01,0.00,0.00,0.00,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.99,0.99,1.00,0.00,0.00,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.97,0.97,0.96,0.96,0.00,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.94,0.94,0.92,0.91,0.90,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.90,0.90,0.89,0.86,0.84,0.82,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.84,0.84,0.83,0.80,0.77,0.75,0.72, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.78,0.78,0.77,0.74,0.72,0.70,0.67, * 0.64,0.00,0.00,0.00,0.00,0.00,0.00, * 0.74,0.74,0.73,0.71,0.68,0.66,0.63, * 0.60,0.59,0.00,0.00,0.00,0.00,0.00, * 0.71,0.71,0.69,0.67,0.64,0.62,0.59, * 0.56,0.55,0.54,0.00,0.00,0.00,0.00, * 0.66,0.66,0.63,0.63,0.60,0.58,0.55, * 0.52,0.50,0.49,0.48,0.00,0.00,0.00, * 0.62,0.62,0.60,0.58,0.56,0.53,0.51, * 0.48,0.46,0.45,0.44,0.43,0.00,0.00, * 0.58,0.58,0.56,0.54,0.52,0.50,0.48, * 0.44,0.43,0.42,0.41,0.40,0.38,0.00, * 0.55,0.55,0.53,0.50,0.48,0.46,0.43, * 0.41,0.40,0.39,0.38,0.36,0.35,0.34/ DATA D5 */1.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 1.00,1.00,0.00,0.00,0.00,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.98,0.98,0.99,0.00,0.00,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.97,0.97,0.96,0.96,0.00,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.94,0.94,0.93,0.92,0.90,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.90,0.90,0.89,0.87,0.84,0.82,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.85,0.85,0.84,0.81,0.78,0.76,0.73, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.79,0.79,0.78,0.75,0.73,0.70,0.68, * 0.65,0.00,0.00,0.00,0.00,0.00,0.00, * 0.75,0.75,0.74,0.72,0.69,0.66,0.64, * 0.61,0.60,0.00,0.00,0.00,0.00,0.00, * 0.72,0.72,0.70,0.68,0.65,0.63,0.60, * 0.57,0.56,0.55,0.00,0.00,0.00,0.00, * 0.68,0.68,0.66,0.64,0.61,0.59,0.56, * 0.53,0.52,0.50,0.49,0.00,0.00,0.00, * 0.64,0.64,0.62,0.60,0.57,0.54,0.52, * 0.49,0.48,0.46,0.45,0.44,0.00,0.00, * 0.61,0.61,0.59,0.56,0.53,0.51,0.49, * 0.46,0.45,0.43,0.42,0.41,0.39,0.00, * 0.58,0.58,0.56,0.53,0.50,0.48,0.46, * 0.43,0.42,0.41,0.39,0.38,0.36,0.36/ DATA D6 */1.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 1.00,1.00,0.00,0.00,0.00,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.98,0.98,0.98,0.00,0.00,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.97,0.96,0.96,0.95,0.00,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.94,0.94,0.93,0.92,0.90,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.90,0.90,0.89,0.87,0.85,0.83,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.86,0.86,0.84,0.82,0.79,0.76,0.74, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.80,0.80,0.79,0.76,0.74,0.71,0.69, * 0.66,0.00,0.00,0.00,0.00,0.00,0.00, * 0.76,0.76,0.75,0.73,0.70,0.67,0.65, * 0.62,0.61,0.00,0.00,0.00,0.00,0.00, * 0.73,0.73,0.72,0.69,0.66,0.64,0.61, * 0.58,0.57,0.56,0.00,0.00,0.00,0.00, * 0.69,0.69,0.68,0.65,0.62,0.60,0.57, * 0.54,0.53,0.52,0.50,0.00,0.00,0.00, * 0.66,0.66,0.64,0.62,0.59,0.56,0.54, * 0.51,0.49,0.48,0.47,0.46,0.00,0.00, * 0.63,0.63,0.61,0.58,0.55,0.53,0.51, * 0.48,0.46,0.45,0.44,0.43,0.41,0.00, * 0.60,0.60,0.58,0.55,0.52,0.50,0.48, * 0.45,0.43,0.42,0.41,0.40,0.39,0.38/ DATA D7 */1.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.99,0.99,0.00,0.00,0.00,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.97,0.97,0.97,0.00,0.00,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.96,0.96,0.95,0.94,0.00,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.94,0.94,0.93,0.91,0.90,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.91,0.91,0.90,0.88,0.85,0.83,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.86,0.86,0.85,0.82,0.80,0.77,0.75, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.81,0.81,0.80,0.77,0.75,0.72,0.70, * 0.67,0.00,0.00,0.00,0.00,0.00,0.00, * 0.77,0.77,0.76,0.74,0.71,0.68,0.66, * 0.63,0.62,0.00,0.00,0.00,0.00,0.00, * 0.74,0.74,0.73,0.70,0.68,0.65,0.62, * 0.60,0.59,0.58,0.00,0.00,0.00,0.00, * 0.71,0.71,0.69,0.67,0.64,0.61,0.59, * 0.56,0.55,0.54,0.52,0.00,0.00,0.00, * 0.68,0.68,0.66,0.64,0.61,0.58,0.56, * 0.53,0.51,0.50,0.48,0.47,0.00,0.00, * 0.65,0.65,0.63,0.60,0.58,0.55,0.52, * 0.50,0.48,0.47,0.45,0.44,0.42,0.00, * 0.63,0.63,0.60,0.57,0.55,0.52,0.49, * 0.47,0.45,0.44,0.43,0.42,0.40,0.39/ DATA D8 */1.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.99,0.99,0.00,0.00,0.00,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.97,0.97,0.97,0.00,0.00,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.96,0.95,0.94,0.93,0.00,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.94,0.94,0.93,0.91,0.90,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.91,0.91,0.90,0.87,0.85,0.83,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.86,0.86,0.85,0.83,0.81,0.78,0.76, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.82,0.82,0.80,0.78,0.75,0.73,0.70, * 0.68,0.00,0.00,0.00,0.00,0.00,0.00, * 0.78,0.78,0.77,0.75,0.72,0.69,0.67, * 0.64,0.63,0.00,0.00,0.00,0.00,0.00, * 0.75,0.75,0.74,0.71,0.69,0.66,0.64, * 0.61,0.60,0.54,0.00,0.00,0.00,0.00, * 0.72,0.72,0.70,0.68,0.65,0.62,0.60, * 0.57,0.56,0.55,0.54,0.00,0.00,0.00, * 0.69,0.69,0.67,0.65,0.62,0.59,0.57, * 0.54,0.52,0.51,0.50,0.48,0.00,0.00, * 0.67,0.67,0.64,0.62,0.59,0.56,0.54, * 0.51,0.50,0.48,0.47,0.46,0.44,0.00, * 0.65,0.65,0.62,0.60,0.57,0.54,0.52, * 0.49,0.48,0.46,0.45,0.44,0.43,0.41/ DATA D9 */1.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.99,0.98,0.00,0.00,0.00,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.97,0.96,0.96,0.00,0.00,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.95,0.94,0.93,0.92,0.00,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.94,0.94,0.92,0.91,0.89,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.91,0.91,0.90,0.87,0.85,0.83,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.87,0.87,0.86,0.84,0.82,0.79,0.77, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.82,0.82,0.81,0.79,0.76,0.74,0.71, * 0.69,0.00,0.00,0.00,0.00,0.00,0.00, * 0.79,0.79,0.78,0.76,0.73,0.71,0.68, * 0.66,0.64,0.00,0.00,0.00,0.00,0.00, * 0.76,0.76,0.75,0.72,0.70,0.67,0.65, * 0.62,0.60,0.59,0.00,0.00,0.00,0.00, * 0.73,0.73,0.72,0.69,0.66,0.63,0.61, * 0.58,0.56,0.55,0.54,0.00,0.00,0.00, * 0.70,0.70,0.68,0.66,0.63,0.60,0.58, * 0.55,0.54,0.52,0.51,0.49,0.00,0.00, * 0.68,0.68,0.66,0.63,0.61,0.58,0.55, * 0.52,0.51,0.49,0.48,0.47,0.45,0.00, * 0.67,0.67,0.64,0.61,0.59,0.56,0.53, * 0.50,0.49,0.47,0.46,0.45,0.44,0.42/ DATA D10 */1.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.98,0.98,0.00,0.00,0.00,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.96,0.96,0.96,0.00,0.00,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.95,0.94,0.93,0.92,0.00,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.94,0.94,0.92,0.91,0.89,0.00,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.91,0.91,0.90,0.88,0.85,0.83,0.00, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.88,0.88,0.87,0.85,0.83,0.80,0.78, * 0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.83,0.83,0.82,0.80,0.77,0.75,0.72, * 0.70,0.00,0.00,0.00,0.00,0.00,0.00, * 0.80,0.80,0.79,0.77,0.74,0.71,0.69, * 0.66,0.64,0.00,0.00,0.00,0.00,0.00, * 0.77,0.77,0.76,0.73,0.71,0.68,0.66, * 0.63,0.62,0.61,0.00,0.00,0.00,0.00, * 0.74,0.74,0.73,0.70,0.67,0.64,0.62, * 0.59,0.58,0.57,0.55,0.00,0.00,0.00, * 0.71,0.70,0.69,0.67,0.64,0.61,0.58, * 0.56,0.54,0.53,0.52,0.50,0.00,0.00, * 0.69,0.69,0.67,0.64,0.62,0.59,0.56, * 0.53,0.52,0.50,0.49,0.48,0.46,0.00, * 0.68,0.68,0.65,0.62,0.60,0.57,0.54, * 0.51,0.50,0.48,0.47,0.46,0.45,0.43/ c********************************************************************* c resolving of problem of zero-values influence on interpolation c********************************************************************* CALL GAP(D1,14,14,1) CALL GAP(D2,14,14,1) CALL GAP(D3,14,14,1) CALL GAP(D4,14,14,1) CALL GAP(D5,14,14,1) CALL GAP(D6,14,14,1) CALL GAP(D7,14,14,1) CALL GAP(D8,14,14,1) CALL GAP(D9,14,14,1) CALL GAP(D10,14,14,1) c********************************************** DO 1 J=1,14 DO 1 I=1,14 F(I,J,1)=D1(I,J) F(I,J,2)=D2(I,J) F(I,J,3)=D3(I,J) F(I,J,4)=D4(I,J) F(I,J,5)=D5(I,J) F(I,J,6)=D6(I,J) F(I,J,7)=D7(I,J) F(I,J,8)=D8(I,J) F(I,J,9)=D9(I,J) F(I,J,10)=D10(I,J) 1 C O N T I N U E RETURN END c********************************************************************* c resolving of problem of zero-values influence on interpolation c smooth extension of Di two-dimension matrix c********************************************************************* C--------------------------------------------------------------- C s m o o t h e x t e n s i o n C (can be used separately for other tasks) C -- 2-dim array C (input - with errors, output - without errors) C , -- array's dimensions C -- averaging radius C -- error's identificator, must be a little higher C--------------------------------------------------------------- SUBROUTINE GAP(X,NX,NY,NPL) DIMENSION X(NX,NY) DATA FUL/0.005/ NA=NPL*2+1 NN=NX IF(NY.GT.NX) NN=NY NN=NN/2+3 IC=NX/2-1 JC=NY/2 1 IC=IC+1 NN=NN+1 IF(X(IC,JC).LT.FUL) GO TO 1 DO 4 K=1,NN II=IC-K JJ=JC-K KK=2*K+1 DO 4 I=1,KK DO 4 L=1,4 IF(L.EQ.1) JR=JJ IF(L.EQ.3) IR=II IF(L.EQ.2) JR=JJ+K*2 IF(L.EQ.4) IR=II+K*2 IF(L.EQ.1.OR.L.EQ.2) IR=II+I-1 IF(L.EQ.3.OR.L.EQ.4) JR=JJ+I-1 IF((IR.LT.1).OR.(IR.GT.NX).OR. * (JR.LT.1).OR.(JR.GT.NY)) GO TO 4 IF(X(IR,JR).GT.FUL) GO TO 4 SU=0. XS=0. DO 2 JA=1,NA DO 2 IA=1,NA ID=IR-NPL+IA-1 JD=JR-NPL+JA-1 IF((ID.LT.1).OR.(ID.GT.NX).OR. * (JD.LT.1).OR.(JD.GT.NY)) GO TO 3 IF(X(ID,JD).LE.FUL) GO TO 3 XS=XS+X(ID,JD) SU=SU+1. 3 CONTINUE 2 CONTINUE IF(SU.LT.0.02) PAUSE ' Error in Interpolation ' IF(SU.LT.0.02) GO TO 4 XS=XS/SU X(IR,JR)=XS 4 CONTINUE DO 5 J=1,NY DO 5 I=1,NX IF(X(I,J).LT.FUL) PAUSE ' Call me ' 5 CONTINUE RETURN END SUBROUTINE WIS (AM,SHR,H) C*************************************************************** C 12.00 (local time) SOLAR ALTITUDE AT 15-TH DAY OF MONTH C--------------------------------------------------------------- C AM - month real input C SHR - latitude real input C H - solar altitude real output C*************************************************************** DIMENSION KD(12) DATA KD/0,31,59,90,120,151,181,212,243,273,304,334/ C*************************************************************** C day from the beginning of year C*************************************************************** J=AM PE=KD(J)+15. C*************************************************************** C solar altitude C*************************************************************** PE1=(PE*0.97826-78.2609)/57.296 H1=(SIN(PE1)*23.75)/57.296 SH=SHR/57.296 H1=sin(sh)*sin(h1)+cos(sh)*cos(h1) H=ASIN(H1)*57.296 IF(H.LE.5.) H=5.01 RETURN END c ------ Dobson rad------------ SUBROUTINE RSWd8 (H,cn,albedo,rdob) dimension a(10),b(10) DATA A/.4,.517,.474,.421,.38,.35,.304,.23,.106,.134/ DATA B/.386,.317,.381,.413,.468,.457,.438,.384,.285,.295/ sol=1368. icn=int(cn) if (h.le.5.01) go to 100 H2=H/57.296 s0=sol*sin(h2) aa=a(icn+1) bb=b(icn+1) tt=aa+bb*sin(h2) rdob=s0*tt*(1.-albedo) go to 200 100 continue rdob=0. 200 continue return end c ------ Reed rad - daily------------ SUBROUTINE RSWr (H,cn,albedo,swri) sol=1353. cc if (cn.lt.9.) rcn=cn+0.06 cc if (cn.ge.9.) rcn=10. rcn=cn/10. if (h.le.5.01) go to 100 H2=H/57.296 s0=sol*sin(h2) r0=s0*(0.61+0.2*(sin(h2))) swri=r0*(1.-0.62*rcn+0.0019*h)*(1.-albedo) go to 200 100 continue swri=0. 200 continue return end C******************************************************************** C s h o r t w a v e r a d i a t i o n C individual values C Sergey Malevsky scheme 1991 -version C Subroutine prepared and tested in January 12-18,1992 St. Petersburg C and in March 13-19,1992 Moscow C by Sergey Malevsky and Sergey Gulev C******************************************************************** C H - solar high deg REAL input C TA - air temperature deg,C REAL input C EZ - air humidity mb REAL " C ANN- low cloudness percent/10 REAL " C ANO- common cloudness percent/10 REAL " C QS - total radiation W/M**2, REAL output C RPOGL-accumulate radiation W/M**2, REAL " C IE - =0, if two levels cloudness INTEGER input C =1, if only common C******************************************************************* SUBROUTINE RSW (H,TA,EZ,ANN,ANO,QS,RPOGL,IE) DIMENSION HSOL(10),ANR1(11),ALB(10,11), *TAR(12),FN(11,11),FN1(11,11,10),FN0(11,10) DATA HSOL/5.,10.,20.,30.,40.,50.,60.,70.,80.,90./ DATA ANR1/0.,1.,2.,3.,4.,5.,6.,7.,8.,9.,10./ DATA TAR/-20.,-15.,-10.,-5.,0.,5.,10.,15.,20.,25.,30.,35./ DATA ALB/0.18,0.20,0.14,0.09,0.06,0.05,0.04,0.04,0.04,0.04, * 0.18,0.20,0.14,0.09,0.06,0.05,0.04,0.04,0.04,0.04, * 0.17,0.19,0.14,0.09,0.06,0.05,0.04,0.04,0.04,0.04, * 0.16,0.18,0.14,0.09,0.06,0.05,0.04,0.04,0.04,0.04, * 0.16,0.18,0.14,0.09,0.06,0.05,0.04,0.04,0.04,0.04, * 0.15,0.17,0.14,0.09,0.06,0.05,0.04,0.04,0.04,0.04, * 0.14,0.15,0.13,0.09,0.07,0.06,0.05,0.05,0.04,0.04, * 0.13,0.14,0.12,0.09,0.07,0.06,0.05,0.05,0.04,0.04, * 0.11,0.12,0.11,0.09,0.07,0.06,0.05,0.05,0.04,0.04, * 0.10,0.10,0.09,0.09,0.07,0.06,0.05,0.05,0.05,0.04, * 0.09,0.09,0.08,0.08,0.07,0.06,0.05,0.05,0.05,0.04/ DATA FN0/1.01,1.02,1.00,0.95,0.95,0.91,0.85,0.78,0.70,0.57,0.37, * 1.00,1.02,1.02,1.00,0.96,0.92,0.85,0.78,0.69,0.56,0.36, * 1.00,1.01,1.01,0.99,0.96,0.92,0.86,0.78,0.68,0.54,0.35, * 1.00,1.00,1.00,0.98,0.96,0.93,0.87,0.80,0.70,0.57,0.38, * 1.00,1.00,0.99,0.98,0.96,0.94,0.88,0.81,0.72,0.60,0.42, * 1.00,0.99,0.98,0.97,0.96,0.94,0.89,0.83,0.75,0.64,0.45, * 1.00,0.98,0.97,0.97,0.96,0.94,0.90,0.84,0.77,0.68,0.49, * 1.00,0.97,0.96,0.96,0.96,0.95,0.92,0.86,0.80,0.71,0.52, * 1.00,0.96,0.96,0.96,0.96,0.96,0.93,0.89,0.83,0.75,0.55, * 1.00,0.96,0.96,0.96,0.96,0.96,0.95,0.91,0.86,0.78,0.59/ IF(H.LE.5.) H=5.01 c********************************************************************* c P2 - constant calculations c 1) P2=f(ez) 2) P2=f(ta) c choose your situation and comment or uncomment others c********************************************************************* c 1) P2=f(ez) - for North Atlantic c********************************************************************* c p2=0.829-0.0078*ez+0.000115*ez*ez c********************************************************************* c 1) P2=f(ez) - for South Atlantic, Indian and Pacific c********************************************************************* ccc p2=0.797-0.0032*ez+0.000034*ez*ez c********************************************************************* c 2) P2=f(ta) - for North Atlantic c********************************************************************* p2=0.799-0.0037*ta c********************************************************************* c 2) P2=f(ta) - for South Atlantic, Indian and Pacific c********************************************************************* ccc p2=0.785-0.0018*ta c********************************************************************* c 2) P2=f(ta) - for World Ocean c********************************************************************* ccc p2=0.790-0.0030*ta c********************************************************************* c C,D - coefficients determination for clear sky incoming radiation c********************************************************************* IF(P2-0.75) 1,2,2 1 C=1.5-2.1*P2+2*P2**2 D=1.808-0.9*P2 GOTO 3 2 C=1.2-1.7*P2+2*P2**2 D=0.38+2.5*P2-2*P2**2 3 CONTINUE H2=H/57.296 IF(D.LT.0.) GO TO 5666 IF(D.EQ.0.) D=0.1 Q0=C*SIN(H2)**D GO TO 78888 5666 Q0=C*1./(SIN(H2)**D) c********************************************************************* 78888 H1=H ifl=0 IF(H1.LE.5.01) ifl=1 c********************************************************************* c albedo calculation -- albedo = F (solar high, total cloud cover) c********************************************************************* CALL INTERP(ANO,H1,HSOL,ANR1,10,11,ALB,ALP) c********************************************************************* c choice the way (TOTAL/TOTAL+LOW) according IE-parameter c********************************************************************* IF(IE.EQ.1) GOTO 16 c********************************************************************* c 1-st way --> T O T A L + L O W c********************************************************************* c preparing (overturning) of 2-levels cloudness matrix (FN1) c********************************************************************* CALL DOF (FN1) c********************************************************************* c cloud cover factor (RF) determination by interpolation c********************************************************************* DO 10 I=1,9 IF(H1.GT.HSOL(I).AND.H1.LE.HSOL(I+1)) GOTO 11 10 CONTINUE 11 N=I DO 12 L=1,11 DO 12 J=1,11 12 FN(L,J)=FN1(J,L,N) CALL INTERP (ANN,ANO,ANR1,ANR1,11,11,FN,S1) N=N+1 DO 13 L=1,11 DO 13 J=1,11 13 FN(L,J)=FN1(J,L,N) CALL INTERP (ANN,ANO,ANR1,ANR1,11,11,FN,S2) RF=S1+(S2-S1)/(HSOL(N)-HSOL(N-1))*(H1-HSOL(N-1)) GOTO 17 c********************************************************************* c 2-nd way --> T O T A L o n l y c********************************************************************* 16 CONTINUE c********************************************************************* c cloud cover factor (RF) determination by interpolation c********************************************************************* CALL INTERP(H1,ANO,ANR1,HSOL,11,10,FN0,RF) 17 CONTINUE c********************************************************************* c incoming shortwave radiation c********************************************************************* QS=Q0*RF*1000. c********************************************************************* c absorbed shortwave radiation c********************************************************************* RPOGL=QS*(1-ALP) c********************************************************************* c small solar high c********************************************************************* if(ifl.eq.0) go to 6789 qs=0. rbl=0. rpogl=0. 6789 continue RETURN END c********************************************************************* c forming of two-level cloud cover coefficints c preparing (overturning) of 2-levels cloudness matrix (F), c including problem of zero-values influence on interpolation c********************************************************************* SUBROUTINE DOF(F) DIMENSION F(11,11,10),D1(11,11),D2(11,11),D3(11,11),D4(11,11), *D5(11,11),D6(11,11),D7(11,11),D8(11,11),D9(11,11),D10(11,11) DATA D1/1.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 1.02,1.02,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 1.01,1.01,1.02,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.98,0.98,0.99,1.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.95,0.95,0.95,0.95,0.95,0.00,0.00,0.00,0.00,0.00,0.00, * 0.91,0.91,0.90,0.90,0.89,0.88,0.00,0.00,0.00,0.00,0.00, * 0.88,0.88,0.87,0.85,0.84,0.82,0.80,0.00,0.00,0.00,0.00, * 0.84,0.84,0.83,0.81,0.78,0.76,0.74,0.72,0.00,0.00,0.00, * 0.80,0.80,0.78,0.76,0.73,0.70,0.67,0.64,0.61,0.00,0.00, * 0.72,0.72,0.70,0.67,0.64,0.61,0.58,0.55,0.52,0.49,0.00, * 0.53,0.53,0.52,0.50,0.48,0.45,0.43,0.41,0.40,0.36,0.34/ DATA D2/1.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 1.01,1.02,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 1.00,1.00,1.02,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.98,0.98,0.99,1.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.96,0.95,0.95,0.95,0.95,0.00,0.00,0.00,0.00,0.00,0.00, * 0.92,0.92,0.91,0.90,0.90,0.89,0.00,0.00,0.00,0.00,0.00, * 0.88,0.88,0.87,0.86,0.84,0.82,0.81,0.00,0.00,0.00,0.00, * 0.84,0.84,0.83,0.80,0.78,0.76,0.74,0.72,0.00,0.00,0.00, * 0.80,0.80,0.78,0.76,0.73,0.70,0.67,0.65,0.61,0.00,0.00, * 0.71,0.71,0.69,0.67,0.64,0.61,0.58,0.56,0.53,0.48,0.00, * 0.52,0.52,0.51,0.48,0.46,0.44,0.42,0.40,0.38,0.36,0.33/ DATA D3/1.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 1.01,1.01,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 1.00,1.00,1.01,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.97,0.97,0.99,0.99,0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.96,0.96,0.95,0.95,0.95,0.00,0.00,0.00,0.00,0.00,0.00, * 0.92,0.92,0.91,0.91,0.91,0.90,0.00,0.00,0.00,0.00,0.00, * 0.89,0.89,0.88,0.87,0.85,0.84,0.83,0.00,0.00,0.00,0.00, * 0.85,0.85,0.84,0.82,0.80,0.77,0.75,0.73,0.00,0.00,0.00, * 0.80,0.80,0.78,0.76,0.73,0.70,0.67,0.64,0.62,0.00,0.00, * 0.71,0.71,0.68,0.66,0.63,0.60,0.57,0.53,0.50,0.50,0.00, * 0.55,0.55,0.52,0.50,0.47,0.44,0.42,0.38,0.37,0.34,0.32/ DATA D4/1.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 1.00,1.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.99,0.99,1.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.97,0.97,0.98,0.98,0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.96,0.96,0.96,0.95,0.95,0.00,0.00,0.00,0.00,0.00,0.00, * 0.93,0.93,0.92,0.92,0.92,0.92,0.00,0.00,0.00,0.00,0.00, * 0.91,0.91,0.90,0.89,0.87,0.86,0.86,0.00,0.00,0.00,0.00, * 0.87,0.87,0.86,0.84,0.82,0.79,0.77,0.75,0.00,0.00,0.00, * 0.82,0.82,0.80,0.78,0.75,0.72,0.69,0.67,0.63,0.00,0.00, * 0.72,0.72,0.71,0.68,0.65,0.62,0.59,0.57,0.54,0.51,0.00, * 0.56,0.56,0.55,0.52,0.50,0.48,0.45,0.43,0.40,0.38,0.35/ DATA D5/1.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 1.00,1.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.98,0.98,0.99,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.97,0.97,0.97,0.98,0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.96,0.96,0.96,0.95,0.95,0.00,0.00,0.00,0.00,0.00,0.00, * 0.94,0.94,0.93,0.93,0.93,0.92,0.00,0.00,0.00,0.00,0.00, * 0.92,0.92,0.91,0.90,0.88,0.87,0.86,0.00,0.00,0.00,0.00, * 0.88,0.88,0.87,0.85,0.83,0.80,0.78,0.76,0.00,0.00,0.00, * 0.84,0.84,0.83,0.80,0.77,0.74,0.71,0.68,0.65,0.00,0.00, * 0.75,0.75,0.74,0.70,0.67,0.64,0.61,0.58,0.55,0.52,0.00, * 0.60,0.60,0.59,0.56,0.53,0.51,0.48,0.45,0.43,0.40,0.38/ DATA D6/1.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.99,0.99,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.98,0.98,0.98,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.97,0.97,0.97,0.97,0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.96,0.96,0.96,0.95,0.95,0.00,0.00,0.00,0.00,0.00,0.00, * 0.95,0.95,0.94,0.93,0.93,0.93,0.00,0.00,0.00,0.00,0.00, * 0.92,0.92,0.92,0.91,0.89,0.88,0.87,0.00,0.00,0.00,0.00, * 0.90,0.90,0.89,0.87,0.84,0.82,0.80,0.78,0.00,0.00,0.00, * 0.86,0.86,0.84,0.82,0.78,0.76,0.73,0.70,0.67,0.00,0.00, * 0.77,0.77,0.76,0.73,0.70,0.67,0.64,0.61,0.58,0.55,0.00, * 0.64,0.64,0.62,0.60,0.57,0.54,0.51,0.49,0.46,0.43,0.40/ DATA D7/1.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.98,0.98,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.97,0.97,0.97,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.96,0.96,0.96,0.97,0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.96,0.96,0.96,0.95,0.95,0.00,0.00,0.00,0.00,0.00,0.00, * 0.95,0.95,0.94,0.94,0.94,0.93,0.00,0.00,0.00,0.00,0.00, * 0.93,0.93,0.92,0.91,0.90,0.88,0.87,0.00,0.00,0.00,0.00, * 0.91,0.91,0.90,0.88,0.86,0.83,0.81,0.79,0.00,0.00,0.00, * 0.87,0.87,0.86,0.83,0.80,0.78,0.75,0.72,0.70,0.00,0.00, * 0.80,0.80,0.78,0.75,0.73,0.69,0.66,0.63,0.60,0.57,0.00, * 0.68,0.68,0.66,0.63,0.60,0.57,0.54,0.51,0.49,0.46,0.43/ DATA D8/1.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.97,0.98,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.96,0.96,0.96,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.97,0.96,0.96,0.96,0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.96,0.96,0.96,0.95,0.95,0.00,0.00,0.00,0.00,0.00,0.00, * 0.95,0.95,0.95,0.94,0.94,0.94,0.00,0.00,0.00,0.00,0.00, * 0.93,0.93,0.93,0.92,0.90,0.88,0.87,0.00,0.00,0.00,0.00, * 0.91,0.91,0.90,0.88,0.86,0.84,0.82,0.80,0.00,0.00,0.00, * 0.89,0.89,0.87,0.85,0.82,0.79,0.77,0.74,0.72,0.00,0.00, * 0.82,0.82,0.81,0.78,0.75,0.72,0.69,0.66,0.63,0.60,0.00, * 0.71,0.71,0.68,0.66,0.63,0.61,0.57,0.54,0.52,0.48,0.45/ DATA D9/1.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.96,0.97,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.96,0.96,0.95,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.95,0.95,0.95,0.96,0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.96,0.96,0.96,0.95,0.95,0.00,0.00,0.00,0.00,0.00,0.00, * 0.96,0.96,0.95,0.94,0.94,0.94,0.00,0.00,0.00,0.00,0.00, * 0.94,0.94,0.93,0.92,0.91,0.90,0.89,0.00,0.00,0.00,0.00, * 0.92,0.92,0.91,0.89,0.87,0.86,0.84,0.82,0.00,0.00,0.00, * 0.90,0.90,0.89,0.86,0.84,0.81,0.78,0.76,0.73,0.00,0.00, * 0.85,0.85,0.83,0.81,0.78,0.75,0.72,0.69,0.66,0.63,0.00, * 0.75,0.75,0.74,0.70,0.67,0.64,0.60,0.57,0.54,0.51,0.48/ DATA D10/1.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.97,0.96,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.96,0.96,0.95,0.00,0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.96,0.96,0.96,0.95,0.00,0.00,0.00,0.00,0.00,0.00,0.00, * 0.96,0.96,0.96,0.95,0.95,0.00,0.00,0.00,0.00,0.00,0.00, * 0.95,0.95,0.95,0.94,0.93,0.91,0.00,0.00,0.00,0.00,0.00, * 0.94,0.94,0.93,0.92,0.91,0.90,0.89,0.00,0.00,0.00,0.00, * 0.93,0.93,0.92,0.90,0.88,0.87,0.85,0.83,0.00,0.00,0.00, * 0.91,0.91,0.90,0.87,0.85,0.82,0.80,0.78,0.75,0.00,0.00, * 0.87,0.87,0.85,0.82,0.80,0.77,0.74,0.71,0.68,0.65,0.00, * 0.78,0.78,0.76,0.73,0.69,0.66,0.62,0.59,0.56,0.52,0.49/ c********************************************************************* c resolving of problem of zero-values influence on interpolation c********************************************************************* CALL GAP(D1,11,11,1) CALL GAP(D2,11,11,1) CALL GAP(D3,11,11,1) CALL GAP(D4,11,11,1) CALL GAP(D5,11,11,1) CALL GAP(D6,11,11,1) CALL GAP(D7,11,11,1) CALL GAP(D8,11,11,1) CALL GAP(D9,11,11,1) CALL GAP(D10,11,11,1) c**************************************************** DO 1 J=1,11 DO 1 I=1,11 F(I,J,1)=D1(I,J) F(I,J,2)=D2(I,J) F(I,J,3)=D3(I,J) F(I,J,4)=D4(I,J) F(I,J,5)=D5(I,J) F(I,J,6)=D6(I,J) F(I,J,7)=D7(I,J) F(I,J,8)=D8(I,J) F(I,J,9)=D9(I,J) F(I,J,10)=D10(I,J) 1 C O N T I N U E RETURN END SUBROUTINE SALT (AD,AM,SHR,DL,SR,IFEWR,H) C*************************************************************** C S O L A R A L T I T U D E C--------------------------------------------------------------- C AD - date real input C AM - month real input C SHR - latitude real input C SHR - longtitude real input C SR - time, GMT real input C IFEWR - number of days real input C in february C H - solar altitude real output C*************************************************************** DIMENSION KD(24) DATA KD/0,31,59,90,120,151,181,212,243,273,304,334,0,31,60, * 91,121,152,182,213,244,274,305,335/ C*************************************************************** C day from the beginning of year C*************************************************************** IF(IFEWR.EQ.28) GOTO 14 J=12+AM GOTO 15 14 J=AM 15 PE=KD(J)+AD C*************************************************************** C r e a l t i m e C*************************************************************** PE2=(PE*1.972+558.028)/57.296 PE1=(PE*0.986+357.014)/57.296 PE1=SIN(PE1)*7.7 PE2=SIN(PE2)*9.5 PE2=(PE1-PE2)*0.0166666 PE1=DL/15.+PE2 IF(DL.LE.0.) T=SR+PE1 IF(DL.GT.0.) T=SR-PE1 197 IF(T.GE.24) GOTO 340 198 IF(T.GE.0.) GOTO 341 T=T+24. GOTO 198 340 T=T-24. GOTO 197 C*************************************************************** C s o l a r a l t i t u d e C*************************************************************** 341 IF(T.GE.12.) GOTO 16 T1=15.0796*(T+12.) GOTO 17 16 T1=15.0796*(T-12.) 17 PE1=(PE*0.97826-78.2609)/57.296 H1=(SIN(PE1)*23.75)/57.296 SH=SHR/57.296 T1 =T1/57.296 H1=sin(SH)*sin(H1)+cos(SH)*cos(H1)*cos(T1) H=ASIN(H1)*57.296 IF(H.LE.5.) H=5.01 RETURN END SUBROUTINE SALT1 (AD,AM,SHR,SR,IFEWR,H) C*************************************************************** C S O L A R A L T I T U D E C--------------------------------------------------------------- C AD - date real input C AM - month real input C SHR - latitude real input C SHR - longtitude real input C SR - time, GMT real input C IFEWR - number of days real input C in february C H - solar altitude real output C*************************************************************** DIMENSION KD(24) DATA KD/0,31,59,90,120,151,181,212,243,273,304,334,0,31,60, * 91,121,152,182,213,244,274,305,335/ C*************************************************************** C day from the beginning of year C*************************************************************** IF(IFEWR.EQ.28) GOTO 14 J=12+AM GOTO 15 14 J=AM 15 PE=KD(J)+AD C*************************************************************** C r e a l t i m e C*************************************************************** c PE2=(PE*1.972+558.028)/57.296 c PE1=(PE*0.986+357.014)/57.296 c PE1=SIN(PE1)*7.7 c PE2=SIN(PE2)*9.5 c PE2=(PE1-PE2)*0.0166666 c PE1=DL/15.+PE2 c IF(DL.LE.0.) T=SR+PE1 c IF(DL.GT.0.) T=SR-PE1 c 197 IF(T.GE.24) GOTO 340 c 198 IF(T.GE.0.) GOTO 341 c T=T+24. c GOTO 198 c 340 T=T-24. c GOTO 197 C*************************************************************** C s o l a r a l t i t u d e C*************************************************************** t=sr 341 IF(T.GE.12.) GOTO 16 T1=15.0796*(T+12.) GOTO 17 16 T1=15.0796*(T-12.) 17 PE1=(PE*0.97826-78.2609)/57.296 H1=(SIN(PE1)*23.75)/57.296 SH=SHR/57.296 T1 =T1/57.296 H1=sin(SH)*sin(H1)+cos(SH)*cos(H1)*cos(T1) H=ASIN(H1)*57.296 ccc IF(H.LE.5.) H=5.01 RETURN END