--- MITgcm_contrib/jscott/igsm/src/radia.F 2006/08/22 20:25:52 1.2 +++ MITgcm_contrib/jscott/igsm/src/radia.F 2007/04/23 21:20:18 1.3 @@ -58,9 +58,10 @@ #if ( defined OCEAN_3D ) #include "AGRID.h" #endif + dimension SWNET(jm0,2),SWIN(jm0,2) #if ( defined CLM ) -#include "CLM.COM" +#include "CLM.h" #endif c @@ -82,9 +83,11 @@ &,BSO4LAND(JM0),BSO4OCEAN(JM0),BSO4TOTAL(JM0) &,bcodmn(JM0,LM0,12) dimension DSWSRF(jm0),DLWSRF(jm0),DSWVIS(jm0),DSWNIR(jm0) + & ,ULWSRF(jm0) integer PCLOUD - common/TSUR/TSURFC(JM0,0:13),TSURFT(JM0),TSURFD(JM0),DTSURF(JM0) - *,cfcld(JM0,3) +! common/TSUR/TSURFC(JM0,0:13),TSURFT(JM0),TSURFD(JM0),DTSURF(JM0) +! *,cfcld(JM0,3) +#include "TSRF.COM" CHARACTER*4 JMNTHF,JMLAST DATA JMLAST /'LAST'/ DATA TF/273.16/,TCIR/258.16/,STBO/.567257E-7/,IFIRST/1/,JDLAST/-9/5045. @@ -174,17 +177,13 @@ print *,' READGHG=',READGHG print *,' CFAEROSOL=',CFAEROSOL #if ( defined PREDICTED_BC) + if(READGHG.eq.0)then +! For runs with black carbon predicted by chemistry model print *,'With black carbon forcing' - print *,' CFBC=',CFBC - nrdbc=0 - read(769),bcodmn -! do nm=1,12 -! do l=1,lm -! do j=1,jm -! bcodmn(j,l,nm)=bcod(1,j,l) -! enddo -! enddo -! enddo + print *,' CFBC=',CFBC + nrdbc=0 + read(769),bcodmn + endif #endif do j=1,jm BSO4TOTAL(j)=1.0 @@ -298,9 +297,9 @@ endif S0X=CFS0X*S0X S0X0=S0X - print *,'Before CALL RADIA0' +! print *,'Before CALL RADIA0' CALL RADIA0 (IO,JM,CO2,READGHG) 5090. - print *,'After CALL RADIA0' +! print *,'After CALL RADIA0' INCHM=NRAD/NDYN 5091. C**** CLOUD LAYER INDICES USED FOR DIAGNOSTICS 5092. DO 43 L=1,LM 5093. @@ -338,12 +337,14 @@ c print *,'From radia S0X0=',S0X0,' S0X=',S0X c print *,'S0X=',S0X,' TNOW=',TNOW ELSE + print *,'From radia TNOW=',TNOW call obssolar(S0X,TNOW) + S0X=S0X/1367. ENDIF ENDIF S0=S0X*1367./RSDIST 5115. S0R=S0 -c print *,'S0=',S0 +! print *,'S0=',S0 C**** CALCULATE AVERAGE COSINE OF ZENITH ANGLE FOR CURRENT COMP3 STEP 5116. C**** AND RADIATION PERIOD 5117. ROT1=TWOPI*TOFDAY/24. 5118. @@ -381,6 +382,7 @@ !ppb(m) to kg per volume base if(READGHG.eq.1) call rtgases(CO2,JMONTH) #if ( defined PREDICTED_BC) + if(READGHG.eq.0)then nrdbc=nrdbc+1 do l=1,lm do j=1,jm @@ -390,6 +392,7 @@ if(nrdbc.eq.12)then nrdbc=0. endif + endif #endif endif JMLAST=JMONTH @@ -411,7 +414,7 @@ DTSURFAV=0. C do j=1,jm - DTSURFAV=DTSURFAV+DTSURF(J)*DXYP(j) + DTSURFAV=DTSURFAV+DT2MGL(J)*DXYP(j) end do !j DTSURFAV=DTSURFAV/AREAG do j=1,jm @@ -428,7 +431,7 @@ print *,'cfcld' print 9456,cfcld print *,' DTSURF' - print 9456,DTSURF + print 9456,DT2MGL print *,' DTSURFAV=',DTSURFAV 9456 format(12f6.2) endif @@ -522,6 +525,10 @@ * JDAY,JYEAR,T,Q,P, * ODATA2,BDATA2,FDATA2,GDATA2,RQT2 & ,CLDSSF,CLDMCF +#if ( defined CLM ) + write(581)asdirclm,asdifclm,aldirclm,aldifclm +#endif + else print *,' NWRCLD=',NWRCLD stop @@ -552,6 +559,10 @@ C endif IF(MODRJ.EQ.0) CALL RCOMPJ 5146. + SWIN(j,1)=0.0 + SWNET(j,1)=0.0 + SWIN(j,2)=0.0 + SWNET(j,2)=0.0 C**** 5147. C**** MAIN I LOOP 5148. C**** 5149. @@ -767,6 +778,10 @@ ENDIF endif ! +#if ( defined IPCC_EMI ) + FULGAS(2)=CO2 +#endif + #if ( defined FIXED_FOR ) FULGAS(2)=1. #endif @@ -854,6 +869,33 @@ C ***** TRSURF(J,ii)=STBO*(POCEAN*TGO**4+POICE*TGOI**4 * +PLICE*TGLI**4+PEARTH*TGE**4)-TRNFLB(1) +! if(TAU.gt.3623.0.and.TAU.lt.3744.0) then + if(ii.eq.-2) then + if(J.eq.29) then +! print *,'From radia TAU=',TAU +! print *,'J=',j,' ii=',ii +! print *,POCEAN,POICE,PLICE,PEARTH +! print *,TGO,TGOI,TGLI,TGE +! print *,' TRNFLB(1)=',TRNFLB(1) +! print *,'TRSURF(J,ii)=',TRSURF(J,ii),' TRDFLB(1)=',TRDFLB(1) +! print *,' TRUFLB(1)=',TRUFLB(1),' SRBO*T**4=', +! & STBO*(POCEAN*TGO**4+POICE*TGOI**4 +! * +PLICE*TGLI**4+PEARTH*TGE**4) +! print *,'CLD',CLDCV,CMC,CSS +! if(ii.eq.2)print *,'TGL=',TGE,STBO*TGE**4 +! write(99,*),'From radia ',TAU,j,ii +! write(99,*), CLDCV,CMC,CSS + write(129), TAU,TOFDAY,CLDCV,TGE,TRSURF(J,ii),SRDFLB(1) +! write(99,*), TAU,TRSURF(J,ii), +! & STBO*(POCEAN*TGO**4+POICE*TGOI**4 +! * +PLICE*TGLI**4+PEARTH*TGE**4), +! & TRNFLB(1) + endif + if(J.eq.30) then + write(130), TAU,TOFDAY,CLDCV,TGE,TRSURF(J,ii),SRDFLB(1) + endif + endif +! endif if(GHSF)then TRSURF(J,ii)=TRSURF(J,ii)+ghostf(1,J) endif @@ -904,6 +946,7 @@ C for TEM CLM DSWSRF(j)=SRDFLB(1) DLWSRF(j)=TRSURF(J,2) + ULWSRF(j)=TRUFLB(1) DSWVIS(j)=SRDVIS DSWNIR(j)=SRDNIR C for TEM CLM @@ -939,9 +982,8 @@ AJ(J,3)=AJ(J,3)+(SRNFLB(1+LM)*COSZ)*POCEAN AJ(J,71)=AJ(J,71)-(TRNFLB(1+LM)-TRNFLB(1))*POCEAN #if ( defined OCEAN_3D ) - solarinc_ocean(J)=solarinc_ocean(J)+SRDFLB(1)*COSZ - solarnet_ocean(J)=solarnet_ocean(J)+SRNFLB(1)*COSZ - navrado(j)=navrado(j)+1 + SWIN(j,1)=SRDFLB(1) + SWNET(j,1)=SRNFLB(1) #endif C DO K=2,9 @@ -962,9 +1004,8 @@ CJ(J,3)=CJ(J,3)+(SRNFLB(1+LM)*COSZ)*POICE CJ(J,71)=CJ(J,71)-(TRNFLB(1+LM)-TRNFLB(1))*POICE #if ( defined OCEAN_3D ) - solarinc_ice(J)=solarinc_ice(J)+SRDFLB(1)*COSZ - solarnet_ice(J)=solarnet_ice(J)+SRNFLB(1)*COSZ - navrad(j)=navrad(j)+1 + SWIN(j,2)=SRDFLB(1) + SWNET(j,2)=SRNFLB(1) #endif C DO K=2,9 @@ -972,40 +1013,40 @@ END DO endif if(CLEAR(J).eq.0)then - PLAND=FDATA(I,J,2) - PWATER=1.-PLAND - POICE=ODATA(I,J,2)*(1.-PLAND) - POCEAN=(1.-PLAND)-POICE - if(POCEAN.LE.1.E-5)then - POCEAN=0. - POICE=PWATER - endif - PLICE=FDATA(I,J,3)*PLAND - PEARTH=PLAND-PLICE - if(ii.eq.1)then - PTYPE=POCEAN - POICE=0. - POCEAN=1. - PLAND=0. - PEARTH=0. - PLICE=0. - else if(ii.eq.3)then - PTYPE=POICE - POICE=1. - POCEAN=0. - PLAND=0. - PEARTH=0. - PLICE=0. - else - PTYPE=PLAND - POCEAN=0. - POICE=0. - PWATER=0. - PLICE=FDATA(I,J,3) - PEARTH=1.-PLICE - PLAND=1. - endif - COSZ=COSZA(I,J) + PLAND=FDATA(I,J,2) + PWATER=1.-PLAND + POICE=ODATA(I,J,2)*(1.-PLAND) + POCEAN=(1.-PLAND)-POICE + if(POCEAN.LE.1.E-5)then + POCEAN=0. + POICE=PWATER + endif + PLICE=FDATA(I,J,3)*PLAND + PEARTH=PLAND-PLICE + if(ii.eq.1)then + PTYPE=POCEAN + POICE=0. + POCEAN=1. + PLAND=0. + PEARTH=0. + PLICE=0. + else if(ii.eq.3)then + PTYPE=POICE + POICE=1. + POCEAN=0. + PLAND=0. + PEARTH=0. + PLICE=0. + else + PTYPE=PLAND + POCEAN=0. + POICE=0. + PWATER=0. + PLICE=FDATA(I,J,3) + PEARTH=1.-PLICE + PLAND=1. + endif + COSZ=COSZA(I,J) if(STRARFOR)then FGOLDU(1)=1.0 elseif(CO2FOR)then @@ -1028,7 +1069,7 @@ endif 568 continue CALL RCOMPX - endif + endif ! CLEAR SRHRCL(J)=SRNFLB(1) TRHRCL(J)=-TRNFLB(1) ALBCL(J)=SRNFLB(1)/(SRDFLB(1)+1.e-20) @@ -1037,6 +1078,12 @@ TRINCL(J)=TRDFLB(1) TRP0CL(J)=TRNFLB(LM+4) TRP1CL(J)=TRNFLB(LM+1) +! if(J.eq.29) then +! if(TAU.gt.3623.0.and.TAU.lt.3744.0) then +! print *,' TRHRCL(j)=',TRHRCL(j) +! print *,' TRINCL(j)=',TRINCL(j) +! endif +! endif COSZ=COSZ2(I,J) do L=LM,1,-1 AJL(J,L,46)=AJL(J,L,46)-TRNFLB(L+1)*PTYPE @@ -1152,15 +1199,34 @@ dlwhr(j)=DLWSRF(j) #endif #if ( defined CLM ) - dsw4clm(j)=DSWSRF(j)*COSZ1(1,j) - dlw4clm(j)=DLWSRF(j) - swinr4clm(j)=DSWNIR(j)*COSZ1(1,j) - swvis4clm(j)=DSWVIS(j)*COSZ1(1,j) + i=1 + dsw4clm(i,j)=DSWSRF(j)*COSZ1(1,j) + dlw4clm(i,j)=DLWSRF(j) + swinr4clm(i,j)=DSWNIR(j)*COSZ1(1,j) + swvis4clm(i,j)=DSWVIS(j)*COSZ1(1,j) +! if(TAU.gt.3633.0.and.TAU.lt.3744.0.and.J.eq.29) then + if(J.eq.-29) then +! print *,TAU,j +! print *,'DLWSRF(j)=',DLWSRF(j), +! & 'dlw4clm(i,j)=',dlw4clm(i,j) + write (229,*)TAU,TOFAY,dsw4clm(i,j),dlw4clm(i,j) + endif + if(J.eq.-30) then + write (230,*)TAU,TOFAY,dsw4clm(i,j),dlw4clm(i,j) + endif c For TEM swtd4tem(j)=swtd4tem(j)+S0*COSZ1(1,j) swsd4tem(j)=swsd4tem(j)+DSWSRF(j)*COSZ1(1,j) nradd4tem(j)=nradd4tem(j)+1 #endif +#if ( defined OCEAN_3D ) + solarinc_ocean(J)=solarinc_ocean(J)+SWIN(j,1)*COSZ1(1,j) + solarnet_ocean(J)=solarnet_ocean(J)+SWNET(j,1)*COSZ1(1,j) + solarinc_ice(J)=solarinc_ice(J)+SWIN(j,2)*COSZ1(1,j) + solarnet_ice(J)=solarnet_ice(J)+SWNET(j,2)*COSZ1(1,j) + navrado(j)=navrado(j)+1 + navrad(j)=navrad(j)+1 +#endif IMAX=IM 5515. IF(J.EQ.1.OR.J.EQ.JM) IMAX=1 5516. if(HPRNT)then