/[MITgcm]/MITgcm_contrib/jscott/igsm/src/radia.F
ViewVC logotype

Diff of /MITgcm_contrib/jscott/igsm/src/radia.F

Parent Directory Parent Directory | Revision Log Revision Log | View Revision Graph Revision Graph | View Patch Patch

revision 1.1 by jscott, Fri Aug 11 19:35:32 2006 UTC revision 1.3 by jscott, Mon Apr 23 21:20:18 2007 UTC
# Line 56  C Line 56  C
56        DIMENSION COE(39)                                                 5043.          DIMENSION COE(39)                                                 5043.  
57        LOGICAL POLE,DC25,HPRNT,WRCLD,CLDFEED                             5044.        LOGICAL POLE,DC25,HPRNT,WRCLD,CLDFEED                             5044.
58  #if ( defined OCEAN_3D )  #if ( defined OCEAN_3D )
59  #include "AGRID.COM"  #include "AGRID.h"
60  #endif  #endif
61           dimension SWNET(jm0,2),SWIN(jm0,2)
62    
63  #if ( defined CLM )  #if ( defined CLM )
64  #include "CLM.COM"  #include "CLM.h"
65  #endif  #endif
66  c  c
67    
# Line 82  c Line 83  c
83       &,BSO4LAND(JM0),BSO4OCEAN(JM0),BSO4TOTAL(JM0)       &,BSO4LAND(JM0),BSO4OCEAN(JM0),BSO4TOTAL(JM0)
84       &,bcodmn(JM0,LM0,12)       &,bcodmn(JM0,LM0,12)
85        dimension DSWSRF(jm0),DLWSRF(jm0),DSWVIS(jm0),DSWNIR(jm0)        dimension DSWSRF(jm0),DLWSRF(jm0),DSWVIS(jm0),DSWNIR(jm0)
86         &   ,ULWSRF(jm0)
87        integer PCLOUD        integer PCLOUD
88        common/TSUR/TSURFC(JM0,0:13),TSURFT(JM0),TSURFD(JM0),DTSURF(JM0)  !     common/TSUR/TSURFC(JM0,0:13),TSURFT(JM0),TSURFD(JM0),DTSURF(JM0)
89       *,cfcld(JM0,3)  !    *,cfcld(JM0,3)
90    #include "TSRF.COM"
91        CHARACTER*4 JMNTHF,JMLAST        CHARACTER*4 JMNTHF,JMLAST
92        DATA JMLAST /'LAST'/        DATA JMLAST /'LAST'/
93        DATA TF/273.16/,TCIR/258.16/,STBO/.567257E-7/,IFIRST/1/,JDLAST/-9/5045.          DATA TF/273.16/,TCIR/258.16/,STBO/.567257E-7/,IFIRST/1/,JDLAST/-9/5045.  
# Line 174  c     print *,' radia.f coefcl=',coefcl Line 177  c     print *,' radia.f coefcl=',coefcl
177        print *,' READGHG=',READGHG        print *,' READGHG=',READGHG
178        print *,' CFAEROSOL=',CFAEROSOL        print *,' CFAEROSOL=',CFAEROSOL
179  #if ( defined PREDICTED_BC)  #if ( defined PREDICTED_BC)
180          if(READGHG.eq.0)then
181    !     For runs with black carbon predicted by chemistry model
182         print *,'With black carbon forcing'         print *,'With black carbon forcing'
183        print *,' CFBC=',CFBC         print *,' CFBC=',CFBC
184        nrdbc=0         nrdbc=0
185        read(769),bcodmn         read(769),bcodmn
186  !     do nm=1,12        endif
 !       do l=1,lm  
 !        do j=1,jm  
 !         bcodmn(j,l,nm)=bcod(1,j,l)  
 !        enddo  
 !       enddo  
 !     enddo  
187  #endif  #endif
188        do j=1,jm        do j=1,jm
189          BSO4TOTAL(j)=1.0          BSO4TOTAL(j)=1.0
# Line 298  c       print *,' 2D ghost forcing' Line 297  c       print *,' 2D ghost forcing'
297        endif        endif
298        S0X=CFS0X*S0X        S0X=CFS0X*S0X
299        S0X0=S0X        S0X0=S0X
300           print *,'Before CALL RADIA0'  !        print *,'Before CALL RADIA0'
301        CALL RADIA0 (IO,JM,CO2,READGHG)                                   5090.          CALL RADIA0 (IO,JM,CO2,READGHG)                                   5090.  
302           print *,'After CALL RADIA0'  !        print *,'After CALL RADIA0'
303           INCHM=NRAD/NDYN                                                5091.             INCHM=NRAD/NDYN                                                5091.  
304  C**** CLOUD LAYER INDICES USED FOR DIAGNOSTICS                          5092.    C**** CLOUD LAYER INDICES USED FOR DIAGNOSTICS                          5092.  
305           DO 43 L=1,LM                                                   5093.             DO 43 L=1,LM                                                   5093.  
# Line 338  C**** RADIATION AVERAGING IN I AND J Line 337  C**** RADIATION AVERAGING IN I AND J
337  c          print *,'From  radia S0X0=',S0X0,' S0X=',S0X  c          print *,'From  radia S0X0=',S0X0,' S0X=',S0X
338  c          print *,'S0X=',S0X,' TNOW=',TNOW  c          print *,'S0X=',S0X,' TNOW=',TNOW
339          ELSE          ELSE
340               print *,'From radia TNOW=',TNOW
341             call obssolar(S0X,TNOW)             call obssolar(S0X,TNOW)
342               S0X=S0X/1367.
343          ENDIF          ENDIF
344        ENDIF        ENDIF
345        S0=S0X*1367./RSDIST                                               5115.          S0=S0X*1367./RSDIST                                               5115.  
346        S0R=S0        S0R=S0
347  c        print *,'S0=',S0  !        print *,'S0=',S0
348  C**** CALCULATE AVERAGE COSINE OF ZENITH ANGLE FOR CURRENT COMP3 STEP   5116.    C**** CALCULATE AVERAGE COSINE OF ZENITH ANGLE FOR CURRENT COMP3 STEP   5116.  
349  C****   AND RADIATION PERIOD                                            5117.    C****   AND RADIATION PERIOD                                            5117.  
350        ROT1=TWOPI*TOFDAY/24.                                             5118.          ROT1=TWOPI*TOFDAY/24.                                             5118.  
# Line 381  c    &           *28.97296245/44.0*1.e-3 Line 382  c    &           *28.97296245/44.0*1.e-3
382                          !ppb(m) to kg per volume base                          !ppb(m) to kg per volume base
383            if(READGHG.eq.1) call rtgases(CO2,JMONTH)            if(READGHG.eq.1) call rtgases(CO2,JMONTH)
384  #if ( defined PREDICTED_BC)  #if ( defined PREDICTED_BC)
385            if(READGHG.eq.0)then
386            nrdbc=nrdbc+1            nrdbc=nrdbc+1
387            do l=1,lm            do l=1,lm
388             do j=1,jm             do j=1,jm
# Line 390  c    &           *28.97296245/44.0*1.e-3 Line 392  c    &           *28.97296245/44.0*1.e-3
392            if(nrdbc.eq.12)then            if(nrdbc.eq.12)then
393              nrdbc=0.              nrdbc=0.
394            endif            endif
395            endif
396  #endif  #endif
397        endif        endif
398        JMLAST=JMONTH        JMLAST=JMONTH
# Line 411  c    &           *28.97296245/44.0*1.e-3 Line 414  c    &           *28.97296245/44.0*1.e-3
414          DTSURFAV=0.          DTSURFAV=0.
415  C  C
416          do j=1,jm          do j=1,jm
417           DTSURFAV=DTSURFAV+DTSURF(J)*DXYP(j)           DTSURFAV=DTSURFAV+DT2MGL(J)*DXYP(j)
418          end do  !j          end do  !j
419           DTSURFAV=DTSURFAV/AREAG           DTSURFAV=DTSURFAV/AREAG
420          do j=1,jm          do j=1,jm
# Line 428  C Line 431  C
431          print *,'cfcld'          print *,'cfcld'
432          print 9456,cfcld          print 9456,cfcld
433          print *,' DTSURF'          print *,' DTSURF'
434          print 9456,DTSURF          print 9456,DT2MGL
435          print *,' DTSURFAV=',DTSURFAV          print *,' DTSURFAV=',DTSURFAV
436   9456   format(12f6.2)   9456   format(12f6.2)
437         endif         endif
# Line 522  CE     END OF READING OF CLOUD Line 525  CE     END OF READING OF CLOUD
525       *              JDAY,JYEAR,T,Q,P,       *              JDAY,JYEAR,T,Q,P,
526       *              ODATA2,BDATA2,FDATA2,GDATA2,RQT2       *              ODATA2,BDATA2,FDATA2,GDATA2,RQT2
527       &             ,CLDSSF,CLDMCF       &             ,CLDSSF,CLDMCF
528    #if ( defined CLM )
529              write(581)asdirclm,asdifclm,aldirclm,aldifclm
530    #endif
531        
532          else          else
533            print *,' NWRCLD=',NWRCLD            print *,' NWRCLD=',NWRCLD
534            stop            stop
# Line 552  C Line 559  C
559  C  C
560        endif        endif
561        IF(MODRJ.EQ.0) CALL RCOMPJ                                        5146.          IF(MODRJ.EQ.0) CALL RCOMPJ                                        5146.  
562          SWIN(j,1)=0.0
563          SWNET(j,1)=0.0
564          SWIN(j,2)=0.0
565          SWNET(j,2)=0.0
566  C****                                                                   5147.    C****                                                                   5147.  
567  C**** MAIN I LOOP                                                       5148.    C**** MAIN I LOOP                                                       5148.  
568  C****                                                                   5149.    C****                                                                   5149.  
# Line 767  c     FGOLDU(3)=XFRADJ*PEARTH Line 778  c     FGOLDU(3)=XFRADJ*PEARTH
778         ENDIF         ENDIF
779        endif        endif
780  !  !
781    #if ( defined IPCC_EMI )
782            FULGAS(2)=CO2
783    #endif
784    
785  #if ( defined FIXED_FOR )  #if ( defined FIXED_FOR )
786          FULGAS(2)=1.          FULGAS(2)=1.
787  #endif  #endif
# Line 854  c Line 869  c
869  C *****  C *****
870        TRSURF(J,ii)=STBO*(POCEAN*TGO**4+POICE*TGOI**4        TRSURF(J,ii)=STBO*(POCEAN*TGO**4+POICE*TGOI**4
871       *  +PLICE*TGLI**4+PEARTH*TGE**4)-TRNFLB(1)       *  +PLICE*TGLI**4+PEARTH*TGE**4)-TRNFLB(1)
872    !     if(TAU.gt.3623.0.and.TAU.lt.3744.0) then
873          if(ii.eq.-2) then
874           if(J.eq.29) then
875    !       print *,'From radia TAU=',TAU
876    !       print *,'J=',j,' ii=',ii
877    !       print *,POCEAN,POICE,PLICE,PEARTH
878    !       print *,TGO,TGOI,TGLI,TGE
879    !       print *,' TRNFLB(1)=',TRNFLB(1)
880    !       print *,'TRSURF(J,ii)=',TRSURF(J,ii),' TRDFLB(1)=',TRDFLB(1)
881    !       print *,' TRUFLB(1)=',TRUFLB(1),' SRBO*T**4=',
882    !    &     STBO*(POCEAN*TGO**4+POICE*TGOI**4
883    !    *     +PLICE*TGLI**4+PEARTH*TGE**4)
884    !        print *,'CLD',CLDCV,CMC,CSS
885    !        if(ii.eq.2)print *,'TGL=',TGE,STBO*TGE**4
886    !        write(99,*),'From radia ',TAU,j,ii
887    !        write(99,*), CLDCV,CMC,CSS
888             write(129), TAU,TOFDAY,CLDCV,TGE,TRSURF(J,ii),SRDFLB(1)
889    !        write(99,*), TAU,TRSURF(J,ii),
890    !    &     STBO*(POCEAN*TGO**4+POICE*TGOI**4
891    !    *     +PLICE*TGLI**4+PEARTH*TGE**4),
892    !    &     TRNFLB(1)
893           endif
894           if(J.eq.30) then
895             write(130), TAU,TOFDAY,CLDCV,TGE,TRSURF(J,ii),SRDFLB(1)
896           endif
897          endif
898    !     endif
899        if(GHSF)then        if(GHSF)then
900           TRSURF(J,ii)=TRSURF(J,ii)+ghostf(1,J)           TRSURF(J,ii)=TRSURF(J,ii)+ghostf(1,J)
901        endif        endif
# Line 904  C ********* Line 946  C *********
946  C        for TEM CLM  C        for TEM CLM
947              DSWSRF(j)=SRDFLB(1)              DSWSRF(j)=SRDFLB(1)
948              DLWSRF(j)=TRSURF(J,2)              DLWSRF(j)=TRSURF(J,2)
949                ULWSRF(j)=TRUFLB(1)
950              DSWVIS(j)=SRDVIS              DSWVIS(j)=SRDVIS
951              DSWNIR(j)=SRDNIR              DSWNIR(j)=SRDNIR
952  C        for TEM CLM  C        for TEM CLM
# Line 939  c         BJ(J,5)=BJ(J,5)+(SRNFLB(1)*COS Line 982  c         BJ(J,5)=BJ(J,5)+(SRNFLB(1)*COS
982            AJ(J,3)=AJ(J,3)+(SRNFLB(1+LM)*COSZ)*POCEAN            AJ(J,3)=AJ(J,3)+(SRNFLB(1+LM)*COSZ)*POCEAN
983            AJ(J,71)=AJ(J,71)-(TRNFLB(1+LM)-TRNFLB(1))*POCEAN            AJ(J,71)=AJ(J,71)-(TRNFLB(1+LM)-TRNFLB(1))*POCEAN
984  #if ( defined OCEAN_3D )  #if ( defined OCEAN_3D )
985            solarinc_ocean(J)=solarinc_ocean(J)+SRDFLB(1)*COSZ            SWIN(j,1)=SRDFLB(1)
986            solarnet_ocean(J)=solarnet_ocean(J)+SRNFLB(1)*COSZ            SWNET(j,1)=SRNFLB(1)
           navrado(j)=navrado(j)+1  
987  #endif  #endif
988  C  C
989            DO  K=2,9            DO  K=2,9
# Line 962  C Line 1004  C
1004            CJ(J,3)=CJ(J,3)+(SRNFLB(1+LM)*COSZ)*POICE            CJ(J,3)=CJ(J,3)+(SRNFLB(1+LM)*COSZ)*POICE
1005            CJ(J,71)=CJ(J,71)-(TRNFLB(1+LM)-TRNFLB(1))*POICE            CJ(J,71)=CJ(J,71)-(TRNFLB(1+LM)-TRNFLB(1))*POICE
1006  #if ( defined OCEAN_3D )  #if ( defined OCEAN_3D )
1007            solarinc_ice(J)=solarinc_ice(J)+SRDFLB(1)*COSZ            SWIN(j,2)=SRDFLB(1)
1008            solarnet_ice(J)=solarnet_ice(J)+SRNFLB(1)*COSZ            SWNET(j,2)=SRNFLB(1)
           navrad(j)=navrad(j)+1  
1009  #endif  #endif
1010  C  C
1011            DO  K=2,9            DO  K=2,9
# Line 972  C Line 1013  C
1013            END DO            END DO
1014            endif            endif
1015        if(CLEAR(J).eq.0)then        if(CLEAR(J).eq.0)then
1016        PLAND=FDATA(I,J,2)                 PLAND=FDATA(I,J,2)        
1017        PWATER=1.-PLAND         PWATER=1.-PLAND
1018        POICE=ODATA(I,J,2)*(1.-PLAND)                     POICE=ODATA(I,J,2)*(1.-PLAND)            
1019        POCEAN=(1.-PLAND)-POICE                       POCEAN=(1.-PLAND)-POICE              
1020        if(POCEAN.LE.1.E-5)then         if(POCEAN.LE.1.E-5)then
1021           POCEAN=0.            POCEAN=0.
1022           POICE=PWATER            POICE=PWATER
1023        endif         endif
1024        PLICE=FDATA(I,J,3)*PLAND                     PLICE=FDATA(I,J,3)*PLAND            
1025        PEARTH=PLAND-PLICE         PEARTH=PLAND-PLICE
1026        if(ii.eq.1)then         if(ii.eq.1)then
1027         PTYPE=POCEAN          PTYPE=POCEAN
1028         POICE=0.          POICE=0.
1029         POCEAN=1.          POCEAN=1.
1030         PLAND=0.          PLAND=0.
1031         PEARTH=0.          PEARTH=0.
1032         PLICE=0.          PLICE=0.
1033        else if(ii.eq.3)then         else if(ii.eq.3)then
1034         PTYPE=POICE          PTYPE=POICE
1035         POICE=1.          POICE=1.
1036         POCEAN=0.          POCEAN=0.
1037         PLAND=0.          PLAND=0.
1038         PEARTH=0.          PEARTH=0.
1039         PLICE=0.          PLICE=0.
1040        else         else
1041         PTYPE=PLAND          PTYPE=PLAND
1042         POCEAN=0.                      POCEAN=0.            
1043         POICE=0.          POICE=0.
1044         PWATER=0.            PWATER=0.  
1045         PLICE=FDATA(I,J,3)          PLICE=FDATA(I,J,3)
1046         PEARTH=1.-PLICE          PEARTH=1.-PLICE
1047         PLAND=1.          PLAND=1.
1048        endif         endif
1049          COSZ=COSZA(I,J)         COSZ=COSZA(I,J)
1050         if(STRARFOR)then         if(STRARFOR)then
1051           FGOLDU(1)=1.0           FGOLDU(1)=1.0
1052         elseif(CO2FOR)then         elseif(CO2FOR)then
# Line 1028  c         print *,'S0 for forcing calcul Line 1069  c         print *,'S0 for forcing calcul
1069         endif         endif
1070    568  continue    568  continue
1071          CALL RCOMPX          CALL RCOMPX
1072        endif        endif  ! CLEAR
1073         SRHRCL(J)=SRNFLB(1)         SRHRCL(J)=SRNFLB(1)
1074         TRHRCL(J)=-TRNFLB(1)         TRHRCL(J)=-TRNFLB(1)
1075         ALBCL(J)=SRNFLB(1)/(SRDFLB(1)+1.e-20)         ALBCL(J)=SRNFLB(1)/(SRDFLB(1)+1.e-20)
# Line 1037  c         print *,'S0 for forcing calcul Line 1078  c         print *,'S0 for forcing calcul
1078         TRINCL(J)=TRDFLB(1)         TRINCL(J)=TRDFLB(1)
1079         TRP0CL(J)=TRNFLB(LM+4)         TRP0CL(J)=TRNFLB(LM+4)
1080         TRP1CL(J)=TRNFLB(LM+1)         TRP1CL(J)=TRNFLB(LM+1)
1081    !     if(J.eq.29) then
1082    !     if(TAU.gt.3623.0.and.TAU.lt.3744.0) then
1083    !       print *,' TRHRCL(j)=',TRHRCL(j)
1084    !       print *,' TRINCL(j)=',TRINCL(j)
1085    !     endif
1086    !     endif
1087         COSZ=COSZ2(I,J)             COSZ=COSZ2(I,J)    
1088         do L=LM,1,-1         do L=LM,1,-1
1089          AJL(J,L,46)=AJL(J,L,46)-TRNFLB(L+1)*PTYPE          AJL(J,L,46)=AJL(J,L,46)-TRNFLB(L+1)*PTYPE
# Line 1152  C**** Line 1199  C****
1199            dlwhr(j)=DLWSRF(j)            dlwhr(j)=DLWSRF(j)
1200  #endif  #endif
1201  #if ( defined CLM )  #if ( defined CLM )
1202            dsw4clm(j)=DSWSRF(j)*COSZ1(1,j)            i=1
1203            dlw4clm(j)=DLWSRF(j)            dsw4clm(i,j)=DSWSRF(j)*COSZ1(1,j)
1204            swinr4clm(j)=DSWNIR(j)*COSZ1(1,j)            dlw4clm(i,j)=DLWSRF(j)
1205            swvis4clm(j)=DSWVIS(j)*COSZ1(1,j)            swinr4clm(i,j)=DSWNIR(j)*COSZ1(1,j)
1206              swvis4clm(i,j)=DSWVIS(j)*COSZ1(1,j)
1207    !         if(TAU.gt.3633.0.and.TAU.lt.3744.0.and.J.eq.29) then
1208              if(J.eq.-29) then
1209    !          print *,TAU,j
1210    !          print *,'DLWSRF(j)=',DLWSRF(j),
1211    !    &     'dlw4clm(i,j)=',dlw4clm(i,j)
1212                write (229,*)TAU,TOFAY,dsw4clm(i,j),dlw4clm(i,j)
1213              endif
1214              if(J.eq.-30) then
1215                write (230,*)TAU,TOFAY,dsw4clm(i,j),dlw4clm(i,j)
1216              endif
1217  c    For TEM  c    For TEM
1218            swtd4tem(j)=swtd4tem(j)+S0*COSZ1(1,j)            swtd4tem(j)=swtd4tem(j)+S0*COSZ1(1,j)
1219            swsd4tem(j)=swsd4tem(j)+DSWSRF(j)*COSZ1(1,j)            swsd4tem(j)=swsd4tem(j)+DSWSRF(j)*COSZ1(1,j)
1220            nradd4tem(j)=nradd4tem(j)+1            nradd4tem(j)=nradd4tem(j)+1
1221  #endif  #endif
1222    #if ( defined OCEAN_3D )
1223              solarinc_ocean(J)=solarinc_ocean(J)+SWIN(j,1)*COSZ1(1,j)
1224              solarnet_ocean(J)=solarnet_ocean(J)+SWNET(j,1)*COSZ1(1,j)
1225              solarinc_ice(J)=solarinc_ice(J)+SWIN(j,2)*COSZ1(1,j)
1226              solarnet_ice(J)=solarnet_ice(J)+SWNET(j,2)*COSZ1(1,j)
1227              navrado(j)=navrado(j)+1
1228              navrad(j)=navrad(j)+1
1229    #endif
1230        IMAX=IM                                                           5515.          IMAX=IM                                                           5515.  
1231        IF(J.EQ.1.OR.J.EQ.JM) IMAX=1                                      5516.          IF(J.EQ.1.OR.J.EQ.JM) IMAX=1                                      5516.  
1232        if(HPRNT)then        if(HPRNT)then
Legend:
Removed from v.1.1  
changed lines
  Added in v.1.3
  ViewVC Help
Powered by ViewVC 1.1.22