!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !! !! !! GNU General Public License !! !! !! !! This file is part of the Flexible Modeling System (FMS). !! !! !! !! FMS is free software; you can redistribute it and/or modify !! !! it and are expected to follow the terms of the GNU General Public !! !! License as published by the Free Software Foundation. !! !! !! !! FMS is distributed in the hope that it will be useful, !! !! but WITHOUT ANY WARRANTY; without even the implied warranty of !! !! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the !! !! GNU General Public License for more details. !! !! !! !! You should have received a copy of the GNU General Public License !! !! along with FMS; if not, write to: !! !! Free Software Foundation, Inc. !! !! 59 Temple Place, Suite 330 !! !! Boston, MA 02111-1307 USA !! !! or see: !! !! http://www.gnu.org/licenses/gpl.txt !! !! !! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! module cloud_zonal_mod !======================================================================= ! ! determines zonal cloud amounts and model levels. ! !======================================================================= use time_manager_mod, only: time_type use time_interp_mod, only: fraction_of_year use fms_mod, only: error_mesg, FATAL, open_file, & close_file, mpp_pe, mpp_root_pe, & write_version_number implicit none private public cloud_zonal, cloud_zonal_init, getcld !------------------- private data used by this module ------------------ integer :: iseason = -1 real, dimension(37,3,4) :: camt real, dimension(37,4) :: phigh,pmidl,ptop,pbtm real :: rad2deg character(len=128) :: version = '$Id: cloud_zonal.F90,v 13.0 2006/03/28 21:07:56 fms Exp $' character(len=128) :: tagname = '$Name: mom4p1_pubrel_dec2009_nnz $' logical :: module_is_initialized = .false. !----------------------------------------------------------------------- !--------cloud amounts every 5 deg. for high,mid, & low----------------- !----------------------------------------------------------------------- !--------------------------- w i n t e r ------------------------------- data camt(:,:,1) & / 0.0850,0.1000,0.1150,0.1300,0.1570,0.1840,0.1970,0.2090, & 0.2030,0.1970,0.1860,0.1740,0.1570,0.1390,0.1360,0.1320,0.1590, & 0.1860,0.2300,0.2300,0.2130,0.1960,0.1760,0.1560,0.1460,0.1360, & 0.1570,0.1780,0.1990,0.2190,0.2420,0.2640,0.2830,0.3010,0.2770, & 0.2520,0.2270,0.0730,0.0780,0.0830,0.0880,0.0990,0.1110,0.1210, & 0.1310,0.1220,0.1120,0.0950,0.0780,0.0690,0.0600,0.0550,0.0500, & 0.0600,0.0700,0.0830,0.0810,0.0760,0.0700,0.0650,0.0590,0.0650, & 0.0700,0.0790,0.0880,0.0990,0.1100,0.1140,0.1180,0.1220,0.1260, & 0.1180,0.1100,0.1020,0.2500,0.2770,0.3040,0.3300,0.3750,0.4200, & 0.4380,0.4550,0.4410,0.4260,0.3840,0.3410,0.2910,0.2410,0.2310, & 0.2200,0.2610,0.3020,0.3500,0.3420,0.3270,0.3120,0.2920,0.2710, & 0.2730,0.2750,0.3340,0.3920,0.4240,0.4560,0.4650,0.4740,0.4830, & 0.4910,0.4710,0.4510,0.4310 / !----------------------------------------------------------------------- !----------------------------------------------------------------------- !--------------------------- s p r i n g ------------------------------- data camt(:,:,2) & / 0.2080,0.2170,0.2260,0.2340,0.2390,0.2440,0.2510,0.2580, & 0.2520,0.2460,0.2370,0.2280,0.2030,0.1780,0.1820,0.1850,0.2100, & 0.2360,0.2630,0.2470,0.2300,0.2120,0.1900,0.1670,0.1770,0.1870, & 0.2020,0.2170,0.2500,0.2830,0.3040,0.3250,0.3330,0.3410,0.2960, & 0.2500,0.2040,0.0720,0.0780,0.0840,0.0890,0.0960,0.1030,0.1180, & 0.1340,0.1270,0.1200,0.1050,0.0900,0.0800,0.0700,0.0680,0.0650, & 0.0690,0.0730,0.0780,0.0760,0.0730,0.0700,0.0670,0.0640,0.0700, & 0.0760,0.0970,0.1180,0.1340,0.1490,0.1460,0.1440,0.1430,0.1420, & 0.1220,0.1020,0.0820,0.3580,0.3750,0.3920,0.4100,0.4080,0.4070, & 0.4050,0.4030,0.3900,0.3770,0.3450,0.3130,0.2760,0.2380,0.2370, & 0.2360,0.2700,0.3050,0.3385,0.3200,0.3030,0.2860,0.2640,0.2410, & 0.2600,0.2780,0.3180,0.3580,0.3970,0.4360,0.4550,0.4740,0.4700, & 0.4650,0.4300,0.3950,0.3600 / !----------------------------------------------------------------------- !----------------------------------------------------------------------- !--------------------------- s u m m e r ------------------------------- data camt(:,:,3) & / 0.2270,0.2520,0.2770,0.3010,0.2830,0.2640,0.2420,0.2190, & 0.1990,0.1780,0.1570,0.1360,0.1460,0.1560,0.1760,0.1960,0.2130, & 0.2300,0.2300,0.1860,0.1590,0.1320,0.1360,0.1390,0.1570,0.1740, & 0.1860,0.1970,0.2030,0.2090,0.1970,0.1840,0.1570,0.1300,0.1150, & 0.1000,0.0850,0.1020,0.1100,0.1180,0.1260,0.1220,0.1180,0.1140, & 0.1100,0.0990,0.0880,0.0790,0.0700,0.0650,0.0590,0.0650,0.0700, & 0.0760,0.0810,0.0830,0.0700,0.0600,0.0500,0.0550,0.0600,0.0690, & 0.0780,0.0950,0.1120,0.1220,0.1310,0.1210,0.1110,0.0990,0.0880, & 0.0830,0.0780,0.0730,0.4310,0.4510,0.4710,0.4910,0.4830,0.4740, & 0.4650,0.4560,0.4240,0.3920,0.3340,0.2750,0.2730,0.2710,0.2920, & 0.3120,0.3270,0.3420,0.3500,0.3020,0.2610,0.2200,0.2310,0.2410, & 0.2910,0.3410,0.3840,0.4260,0.4410,0.4550,0.4380,0.4200,0.3750, & 0.3300,0.3040,0.2770,0.2500 / !----------------------------------------------------------------------- !----------------------------------------------------------------------- !----------------------------- f a l l --------------------------------- data camt(:,:,4) & / 0.2040,0.2500,0.2960,0.3410,0.3330,0.3250,0.3040,0.2830, & 0.2500,0.2170,0.2020,0.1870,0.1770,0.1670,0.1900,0.2120,0.2300, & 0.2470,0.2630,0.2360,0.2100,0.1850,0.1820,0.1780,0.2030,0.2280, & 0.2370,0.2460,0.2520,0.2580,0.2510,0.2440,0.2390,0.2340,0.2260, & 0.2170,0.2080,0.0820,0.1020,0.1220,0.1420,0.1430,0.1440,0.1460, & 0.1490,0.1340,0.1180,0.0970,0.0760,0.0700,0.0640,0.0670,0.0700, & 0.0730,0.0760,0.0780,0.0730,0.0690,0.0650,0.0680,0.0700,0.0800, & 0.0900,0.1050,0.1200,0.1270,0.1340,0.1180,0.1030,0.0960,0.0890, & 0.0840,0.0780,0.0720,0.3600,0.3950,0.4300,0.4650,0.4700,0.4740, & 0.4550,0.4360,0.3970,0.3580,0.3180,0.2780,0.2600,0.2410,0.2640, & 0.2860,0.3030,0.3200,0.3385,0.3050,0.2700,0.2360,0.2370,0.2380, & 0.2760,0.3130,0.3450,0.3770,0.3900,0.4030,0.4050,0.4070,0.4080, & 0.4100,0.3920,0.3750,0.3580 / !----------------------------------------------------------------------- !----------------------------------------------------------------------- !---------pressure for high clouds (winter,spring,summer,fall)---------- data phigh & /42253., 41101., 39434., 37223., 35249., 34212., 33924., 33682., & 32911., 31304., 28946., 26365., 23980., 21839., 20149., 19267., & 19074., 19104., 19115., 19111., 19113., 19126., 19191., 19445., & 20168., 21652., 23863., 26357., 28787., 31081., 32924., 33884., & 34093., 34061., 34049., 34053., 34053., & 34053., 34053., 34053., 34051., 34038., 33967., 33690., 32907., & 31304., 28947., 26364., 23980., 21839., 20149., 19267., 19074., & 19104., 19115., 19111., 19113., 19126., 19191., 19445., 20168., & 21652., 23863., 26357., 28787., 31081., 32924., 33884., 34093., & 34061., 34049., 34053., 34053., 34053., & 34053., 34053., 34051., 34038., 33967., 33690., 32907., 31304., & 28947., 26364., 23980., 21839., 20149., 19267., 19074., 19104., & 19115., 19111., 19113., 19126., 19191., 19445., 20168., 21652., & 23863., 26357., 28787., 31081., 32926., 33900., 34184., 34442., & 35242., 36864., 39075., 41055., 42253., & 34053., 34053., 34053., 34051., 34038., 33967., 33690., 32907., & 31304., 28947., 26364., 23980., 21839., 20149., 19267., 19074., & 19104., 19115., 19111., 19113., 19126., 19191., 19445., 20168., & 21652., 23863., 26357., 28787., 31081., 32924., 33884., 34093., & 34061., 34049., 34053., 34053., 34053./ !----------------------------------------------------------------------- !----------------------------------------------------------------------- !--------pressure for middle clouds (winter,spring,summer,fall)--------- data pmidl & /63172., 63172., 63172., 63172., 63172., 63172., 63172., 63172., & 63172., 63172., 63172., 63172., 63172., 63170., 63155., 63080., & 62786., 61964., 60324., 58078., 56082., 55059., 54916., 55244., & 56080., 57715., 59961., 61959., 62996., 63215., 63180., 63167., & 63172., 63172., 63172., 63172., 63172., & 63172., 63172., 63172., 63172., 63172., 63172., 63172., 63170., & 63155., 63080., 62786., 61964., 60324., 58078., 56079., 55043., & 54824., 54859., 54874., 54883., 54959., 55252., 56075., 57715., & 59961., 61959., 62996., 63215., 63180., 63167., 63172., 63172., & 63172., 63172., 63172., 63172., 63172., & 63172., 63172., 63172., 63172., 63172., 63170., 63155., 63080., & 62786., 61964., 60324., 58078., 56082., 55059., 54916., 55244., & 56080., 57715., 59961., 61959., 62996., 63215., 63180., 63167., & 63172., 63172., 63172., 63172., 63172., 63172., 63172., 63172., & 63172., 63172., 63172., 63172., 63172., & 63172., 63172., 63172., 63172., 63172., 63172., 63172., 63170., & 63155., 63080., 62786., 61964., 60324., 58078., 56079., 55043., & 54824., 54859., 54874., 54883., 54959., 55252., 56075., 57715., & 59961., 61959., 62996., 63215., 63180., 63167., 63172., 63172., & 63172., 63172., 63172., 63172., 63172./ !----------------------------------------------------------------------- !----------------------------------------------------------------------- !-------pressure for top of low clouds (winter,spring,summer,fall)------ data ptop & /75056., 75056., 75056., 75056., 75056., 75056., 75055., 75054., & 75040., 74969., 74694., 73923., 72386., 70281., 68406., 67421., & 67141., 66886., 66099., 64550., 62622., 61462., 62058., 64061., & 66082., 67189., 67675., 68412., 69937., 72047., 73919., 74890., & 75096., 75063., 75051., 75056., 75056., & 75056., 75056., 75056., 75055., 75054., 75040., 74969., 74694., & 73923., 72386., 70281., 68408., 67437., 67232., 67264., 67276., & 67272., 67272., 67272., 67272., 67272., 67274., 67287., 67358., & 67633., 68404., 69941., 72046., 73919., 74890., 75096., 75063., & 75051., 75056., 75056., 75055., 75056., & 75056., 75055., 75054., 75040., 74969., 74694., 73923., 72386., & 70279., 68392., 67347., 66853., 66087., 64555., 62622., 61461., & 62058., 64059., 66068., 67118., 67400., 67641., 68400., 69942., & 72047., 73919., 74890., 75096., 75063., 75051., 75056., 75056., & 75055., 75056., 75056., 75056., 75056., & 75056., 75056., 75056., 75055., 75054., 75040., 74969., 74694., & 73923., 72386., 70281., 68408., 67437., 67232., 67264., 67276., & 67272., 67272., 67272., 67272., 67272., 67274., 67287., 67358., & 67633., 68404., 69941., 72046., 73919., 74890., 75096., 75063., & 75051., 75056., 75056., 75055., 75056./ !----------------------------------------------------------------------- !----------------------------------------------------------------------- !------pressure for bottom of low clouds (winter,spring,summer,fall)---- data pbtm & /93819., 93818., 93809., 93761., 93574., 93051., 92008., 90580., & 89309., 88650., 88510., 88532., 88541., 88538., 88538., 88538., & 88538., 88538., 88538., 88538., 88538., 88538., 88538., 88538., & 88538., 88538., 88538., 88538., 88539., 88548., 88596., 88783., & 89306., 90348., 91772., 93048., 93819., & 93819., 93575., 93050., 92008., 90580., 89309., 88650., 88510., & 88532., 88541., 88538., 88538., 88538., 88538., 88538., 88538., & 88538., 88538., 88538., 88538., 88538., 88538., 88538., 88538., & 88538., 88538., 88538., 88539., 88548., 88596., 88783., 89306., & 90349., 91778., 93048., 93699., 93819., & 93819., 93077., 92003., 90579., 89309., 88650., 88510., 88532., & 88541., 88538., 88538., 88538., 88538., 88538., 88538., 88538., & 88538., 88538., 88538., 88538., 88538., 88538., 88538., 88538., & 88539., 88548., 88596., 88783., 89306., 90349., 91777., 93048., & 93707., 93847., 93825., 93817., 93819., & 93819., 93575., 93050., 92008., 90580., 89309., 88650., 88510., & 88532., 88541., 88538., 88538., 88538., 88538., 88538., 88538., & 88538., 88538., 88538., 88538., 88538., 88538., 88538., 88538., & 88538., 88538., 88538., 88539., 88548., 88596., 88783., 89306., & 90349., 91778., 93048., 93699., 93819./ !----------------------------------------------------------------------- ! note: the fels-schwarzkopf radiation code permits bi-spectral ! cloud reflectivity associated with cloud cdwtr droplets: ! --> visible band - (cao3sw and cuvrf); ! --> near infra-red band - (cah2sw and cirrf). ! the f-s code assumes that all gaseous absorption by ! cdwtr vapor occurs in the near infra-red band. ! thus, original code contains cbsw and cirab. ! we shall include cbo3sw and cuvab and let cbsw = cbh2sw. ! however, these spectral absorptivities will be set to zero. real, dimension(3) :: cao3sw = (/ 0.210, 0.450, 0.590 /) real, dimension(3) :: cah2sw = (/ 0.210, 0.450, 0.590 /) real, dimension(3) :: cbsw = (/ 0.005, 0.020, 0.035 /) !----------------------------------------------------------------------- contains !####################################################################### subroutine cloud_zonal_init (season) !----------------------------------------------------------------------- ! ! initialization routine for retrieval of ! zonal cloud amounts and level indices. ! ! input argument ! -------------- ! ! season scalar integer between 1-5 ! where 1-4 uses fixed data (1=winter, 2=spring, etc.) ! season=5 is seasonal varying clouds ! integer, intent(in) :: season character(len=34) :: err_string !----------------------------------------------------------------------- !---- error checks ----- if (iseason /= -1) then return !DEL call error_mesg ('cloud_zonal_init', & !DEL 'initialization routine can not be called twice.', FATAL) endif if (season < 1 .or. season > 5) then write (err_string,9001) season 9001 format ('invalid value of season=',i10) call error_mesg ('cloud_zonal_init', err_string, FATAL) endif iseason=season rad2deg=90./acos(0.0) !---- print version number to logfile ---- if ( mpp_pe() == mpp_root_pe() ) then call write_version_number(version, tagname) endif module_is_initialized = .true. !----------------------------------------------------------------------- end subroutine cloud_zonal_init !####################################################################### subroutine cloud_zonal_end module_is_initialized =.false. !----------------------------------------------------------------------- end subroutine cloud_zonal_end !####################################################################### subroutine getcld (time, lat, phalf, ktopsw, kbtmsw, cldamt) !----------------------------------------------------------------------- ! ! routine for retrieval of zonal cloud amounts and level indices. ! ! input arguments ! -------------- ! ! time time of year (time_type) ! lat latitudes in radians, dimensioned by ncol ! phalf pressure at model layer interfaces, ! dimensioned mxcol x nlev, although only ! the first ncol points of the first dimension ! are processed ! ! output arguments ! ---------------- ! ! (all output arguments are dimensioned ncol x 3; the second ! dimension represents high, middle, and low clouds) ! ! ktopsw model layer interface indices for cloud tops ! kbtmsw model layer interface indices for cloud bottoms ! cldamt fractional cloud amounts ! type(time_type), intent(in) :: time real, intent(in) :: lat(:,:), phalf(:,:,:) integer, intent(out) :: ktopsw(:,:,:),kbtmsw(:,:,:) real , intent(out), optional :: cldamt(:,:,:) !----------------------------------------------------------------------- real dlag,rsea,dsea,than,pptop,ppmid,ppbtm,phdeg, & pcnt,ppup,pplo, fyear integer i,j,k,n,isea1,isea2,kmx,lat1,ltop,lbtm !----------------------------------------------------------------------- real cam (37,3),phi(37),pmi(37),plu(37),plb(37) !----------------------------------------------------------------------- if (iseason == -1) then call cloud_zonal_init (5) endif !----------------------------------------------------------------------- !--------------time interpolation if seasonally varying----------------- !----------------------------------------------------------------------- if (iseason == 5) then fyear = fraction_of_year (time) dlag=1./24. rsea=4.*(fyear-dlag)+1.0 isea1=int(rsea) dsea=rsea-float(isea1) if (isea1 < 1) isea1=isea1+4 if (isea1 > 4) isea1=isea1-4 isea2=isea1+1 if (isea2 > 4) isea2=isea2-4 cam(:,:) = camt(:,:,isea1)+dsea*(camt(:,:,isea2)-camt(:,:,isea1)) phi(:) = phigh(:,isea1) + dsea*(phigh(:,isea2)-phigh(:,isea1)) pmi(:) = pmidl(:,isea1) + dsea*(pmidl(:,isea2)-pmidl(:,isea1)) plu(:) = ptop (:,isea1) + dsea*(ptop (:,isea2)-ptop (:,isea1)) plb(:) = pbtm (:,isea1) + dsea*(pbtm (:,isea2)-pbtm (:,isea1)) else !----------------------------------------------------------------------- !--------------no interpolation: not seasonally varying----------------- !----------------------------------------------------------------------- cam(:,:)=camt(:,:,iseason) phi(:) = phigh(:,iseason) pmi(:) = pmidl(:,iseason) plu(:) = ptop (:,iseason) plb(:) = pbtm (:,iseason) endif !----------------------------------------------------------------------- !======================================================================= !---------------------vertical interpolation---------------------------- !----------------------------------------------------------------------- !----------------------------------------------------------------------- ! cloud interpolation ---> 18 level data to more than 18 levels ! ! high cloud (lc=1) interpolates to upper layer interfaces ! mid cloud (lc=2) interpolates to middle of layer ! low cloud (lc=3) interpolates to upper and lower ! layer interfaces ! ! **** note: high and middle clouds are one layer thick !----------------------------------------------------------------------- kmx=size(phalf,3) do n=1,3 do j=1,size(lat,2) do i=1,size(lat,1) !----------------------------------------------------------------------- phdeg=max(min(lat(i,j)*rad2deg,90.),-90.) lat1=19.000-0.20*phdeg lat1=max(min(lat1,36),1) than=(19-lat1)-0.20*phdeg !----------------------------------------------------------------------- ! ------ cloud amount ------ if (present(cldamt)) then cldamt(i,j,n)=cam(lat1,n)+than*(cam(lat1+1,n)-cam(lat1,n)) endif !******************** high and low clouds ****************************** if (n == 1 .or. n == 3) then if (n == 1) pptop=phi(lat1)+than*(phi(lat1+1)-phi(lat1)) if (n == 3) pptop=plu(lat1)+than*(plu(lat1+1)-plu(lat1)) ! --- upper interface --- do 801 k=2,kmx if (pptop > phalf(i,j,k)) go to 801 if ( abs(pptop-phalf(i,j,k-1)) < abs(pptop-phalf(i,j,k)) ) then ltop=k-1 else ltop=k endif go to 802 801 continue ltop=kmx 802 continue if (n == 1) lbtm=ltop+1 endif !*********************************************************************** ! !******************** middle clouds only ******************************* if (n == 2) then ppmid=pmi(lat1)+than*(pmi(lat1+1)-pmi(lat1)) ! --- upper interface --- do 803 k=2,kmx if (ppmid > phalf(i,j,k)) go to 803 ltop=k-1 lbtm=k go to 804 803 continue ltop=kmx lbtm=kmx 804 continue endif !*********************************************************************** ! !******************** low clouds only********************************** if (n == 3) then ppbtm=plb(lat1)+than*(plb(lat1+1)-plb(lat1)) ! --- lower interface --- do 807 k=2,kmx if (ppbtm > phalf(i,j,k)) go to 807 if ( abs(ppbtm-phalf(i,j,k-1)) < abs(ppbtm-phalf(i,j,k)) ) then lbtm=k-1 else lbtm=k endif go to 808 807 continue lbtm=kmx 808 continue if (ltop == lbtm) then pcnt=pptop+ppbtm ppup=phalf(i,j,ltop-1)+phalf(i,j,lbtm ) pplo=phalf(i,j,ltop )+phalf(i,j,lbtm+1) if ( abs(pcnt-ppup) <= abs(pcnt-pplo) ) then ltop=ltop-1 else lbtm=lbtm+1 endif endif endif !*********************************************************************** if ( n >= 2 ) then if ( ltop < kbtmsw(i,j,n-1) ) ltop=kbtmsw(i,j,n-1) if ( lbtm <= ltop ) lbtm=ltop+1 endif ktopsw(i,j,n)=ltop kbtmsw(i,j,n)=lbtm !----------------------------------------------------------------------- enddo enddo enddo !----------------------------------------------------------------------- end subroutine getcld !####################################################################### subroutine cloud_zonal (time, lat, phalf, & nclds, ktopsw, kbtmsw, ktoplw, kbtmlw, & cldamt, cuvrf, cirrf, cirab, emcld) !----------------------------------------------------------------------- type(time_type), intent(in) :: time real, intent(in) :: lat(:,:), phalf(:,:,:) integer, intent(out), dimension(:,:) :: nclds integer, intent(out), dimension(:,:,:) :: ktopsw,kbtmsw,ktoplw,kbtmlw real, intent(out), dimension(:,:,:) :: cldamt,cuvrf,cirrf,cirab,emcld !----------------------------------------------------------------------- integer,dimension(size(ktopsw,1),size(ktopsw,2),3) :: ktopsw3,kbtmsw3 real,dimension(size(cldamt,1),size(cldamt,2),3) :: cldamt3 integer k !----------------------------------------------------------------------- !! kp1=size(ktopsw,3) ! ------- default clouds -------- !! cldamt=0.0; emcld =1.0 !! cuvrf =0.0; cirrf =0.0; cirab =0.0 !! ktopsw=kp1; kbtmsw=kp1 !! ktoplw=kp1; kbtmlw=kp1 ! ---- reset top properties ---- !! ktopsw(:,:,1)=1 !! kbtmsw(:,:,1)=0 !! ktoplw(:,:,1)=1 !! kbtmlw(:,:,1)=0 !-------- insert the 3 cloud layers into array positions 1-3 ----------- call getcld (time, lat, phalf, ktopsw3, kbtmsw3, cldamt3) nclds = 3 emcld = 1.0 do k = 1,3 ktopsw(:,:,k)=ktopsw3(:,:,k) kbtmsw(:,:,k)=kbtmsw3(:,:,k) cldamt(:,:,k)=cldamt3(:,:,k) ktoplw(:,:,k)=ktopsw(:,:,k) kbtmlw(:,:,k)=kbtmsw(:,:,k)-1 cuvrf (:,:,k)=cao3sw(k) cirrf (:,:,k)=cah2sw(k) cirab (:,:,k)=cbsw (k) enddo !----------------------------------------------------------------------- end subroutine cloud_zonal !####################################################################### end module cloud_zonal_mod