!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !! !! !! 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 LONGWAVE_MOD !----------------------------------------------------------------------- USE RDPARM_MOD, ONLY: LMAX USE RDPARM_MOD, ONLY: LM1,LP1,LP2,LL,LLP1,LLM1,LP1M,LP1V,LL3P USE RDPARM_MOD, ONLY: NBLW,NBLX,NBLY,NBLM,INLTE,INLTEP,NNLTE USE HCONST_MOD, ONLY: DIFFCTR,GINV,P0,P0INV,GP0INV,P0XZP2,P0XZP8 USE HCONST_MOD, ONLY: RADCON,RADCON1,RATH2OMW,SECPDA Use FMS_Mod, ONLY: Error_Mesg, FATAL, NOTE, mpp_pe, & mpp_root_pe, write_version_number Use CO2_Data_Mod, ONLY: CO251,CO258,CDT51,CDT58,C2D51,C2D58, & CO2M51,CO2M58,CDTM51,CDTM58,C2DM51, & C2DM58, STEMP,GTEMP, B0,B1,B2,B3 Use CO2_Data_Mod, ONLY: CO231,CO238,CDT31,CDT38,C2D31,C2D38 Use CO2_Data_Mod, ONLY: CO271,CO278,CDT71,CDT78,C2D71,C2D78 Use CO2_Data_Mod, ONLY: CO211,CO218 ! ----------------------------------------------------------- implicit none private !----------------------------------------------------------------------- !--------------------- G L O B A L D A T A --------------------------- !----------------------------------------------------------------------- ! ! Random band parameters for the longwave calcualtions using ! 10 cm-1 wide bands. The 15 um co2 complex is 2 bands, ! 560-670 and 670-800 cm-1. Ozone coefficients are in 3 bands, ! 670-800 (14.1 um), 990-1070 and 1070-1200 (9.6 um). ! The (NBLW) bands now include: ! ! 56 BANDS, 10 CM-1 WIDE 0 - 560 CM-1 ! 2 BANDS, 15 UM COMPLEX 560 - 670 CM-1 ! 670 - 800 CM-1 ! 3 "CONTINUUM" BANDS 800 - 900 CM-1 ! 900 - 990 CM-1 ! 1070 - 1200 CM-1 ! 1 BAND FOR 9.6 UM BAND 990 - 1070 CM-1 ! 100 BANDS, 10 CM-1 WIDE 1200 - 2200 CM-1 ! 1 BAND FOR 4.3 UM SRC 2270 - 2380 CM-1 ! ! Thus NBLW presently equals 163 ! All bands are arranged in order of increasing wavenumbers. ! ! ARNDM = RANDOM "A" PARAMETER FOR (NBLW) BANDS ! BRNDM = RANDOM "B" PARAMETER FOR (NBLW) BANDS ! BETAD = CONTINUUM COEFFICIENTS FOR (NBLW) BANDS ! AP,BP = CAPPHI COEFFICIENTS FOR (NBLW) BANDS ! ATP,BTP = CAPPSI COEFFICIENTS FOR (NBLW) BANDS ! BANDLO = LOWEST FREQUENCY IN EACH OF (NBLW) FREQ. BANDS ! BANDHI = HIGHEST FREQUENCY IN EACH OF (NBLW) FREQ. BANDS ! AO3RND = RANDOM "A" PARAMETER FOR OZONE IN (3) OZONE ! BANDS. ! BO3RND = RANDOM "B" PARAMETER FOR OZONE IN (3) OZONE ! BANDS ! AB15 = THE PRODUCT ARNDM*BRNDM FOR THE TWO BANDS ! REPRESENTING THE 15 UM BAND COMPLEX OF CO2 ! ! Data for ARNDM,BRNDM,AP,BP,ATP,BTP,AO3RND,BO3RND are obtained ! by using the AFGL 1982 catalog. Continuum coefficients are from ! Roberts (1976). This data was formerly in COMMON /BANDTA/. REAL ARNDM(NBLW),BRNDM(NBLW),BETAD(NBLW),AP(NBLW), & BP(NBLW),ATP(NBLW),BTP(NBLW),BANDLO(NBLW), & BANDHI(NBLW),AO3RND(3),BO3RND(3),AB15(2) ! ---------------------------------------------------------- ! The following data statements are band parameters obtained ! using the 1982 AFGL catalog on the specified bands. ! ---------------------------------------------------------- DATA AO3RND /0.543368E+02, 0.234676E+04, 0.384881E+02/ DATA BO3RND /0.526064E+01, 0.922424E+01, 0.496515E+01/ integer :: i DATA (ARNDM(i),i=1,64) / & 0.354693E+00, 0.269857E+03, 0.167062E+03, 0.201314E+04, & 0.964533E+03, 0.547971E+04, 0.152933E+04, 0.599429E+04, & 0.699329E+04, 0.856721E+04, 0.962489E+04, 0.233348E+04, & 0.127091E+05, 0.104383E+05, 0.504249E+04, 0.181227E+05, & 0.856480E+03, 0.136354E+05, 0.288635E+04, 0.170200E+04, & 0.209761E+05, 0.126797E+04, 0.110096E+05, 0.336436E+03, & 0.491663E+04, 0.863701E+04, 0.540389E+03, 0.439786E+04, & 0.347836E+04, 0.130557E+03, 0.465332E+04, 0.253086E+03, & 0.257387E+04, 0.488041E+03, 0.892991E+03, 0.117148E+04, & 0.125880E+03, 0.458852E+03, 0.142975E+03, 0.446355E+03, & 0.302887E+02, 0.394451E+03, 0.438112E+02, 0.348811E+02, & 0.615503E+02, 0.143165E+03, 0.103958E+02, 0.725108E+02, & 0.316628E+02, 0.946456E+01, 0.542675E+02, 0.351557E+02, & 0.301797E+02, 0.381010E+01, 0.126319E+02, 0.548010E+01, & 0.600199E+01, 0.640803E+00, 0.501549E-01, 0.167961E-01, & 0.178110E-01, 0.170166E+00, 0.273514E-01, 0.983767E+00/ DATA (ARNDM(i),i=65,128) / & 0.753946E+00, 0.941763E-01, 0.970547E+00, 0.268862E+00, & 0.564373E+01, 0.389794E+01, 0.310955E+01, 0.128235E+01, & 0.196414E+01, 0.247113E+02, 0.593435E+01, 0.377552E+02, & 0.305173E+02, 0.852479E+01, 0.116780E+03, 0.101490E+03, & 0.138939E+03, 0.324228E+03, 0.683729E+02, 0.471304E+03, & 0.159684E+03, 0.427101E+03, 0.114716E+03, 0.106190E+04, & 0.294607E+03, 0.762948E+03, 0.333199E+03, 0.830645E+03, & 0.162512E+04, 0.525676E+03, 0.137739E+04, 0.136252E+04, & 0.147164E+04, 0.187196E+04, 0.131118E+04, 0.103975E+04, & 0.621637E+01, 0.399459E+02, 0.950648E+02, 0.943161E+03, & 0.526821E+03, 0.104150E+04, 0.905610E+03, 0.228142E+04, & 0.806270E+03, 0.691845E+03, 0.155237E+04, 0.192241E+04, & 0.991871E+03, 0.123907E+04, 0.457289E+02, 0.146146E+04, & 0.319382E+03, 0.436074E+03, 0.374214E+03, 0.778217E+03, & 0.140227E+03, 0.562540E+03, 0.682685E+02, 0.820292E+02, & 0.178779E+03, 0.186150E+03, 0.383864E+03, 0.567416E+01/ DATA (ARNDM(i),i=129,163) / & 0.225129E+03, 0.473099E+01, 0.753149E+02, 0.233689E+02, & 0.339802E+02, 0.108855E+03, 0.380016E+02, 0.151039E+01, & 0.660346E+02, 0.370165E+01, 0.234169E+02, 0.440206E+00, & 0.615283E+01, 0.304077E+02, 0.117769E+01, 0.125248E+02, & 0.142652E+01, 0.241831E+00, 0.483721E+01, 0.226357E-01, & 0.549835E+01, 0.597067E+00, 0.404553E+00, 0.143584E+01, & 0.294291E+00, 0.466273E+00, 0.156048E+00, 0.656185E+00, & 0.172727E+00, 0.118349E+00, 0.141598E+00, 0.588581E-01, & 0.919409E-01, 0.155521E-01, 0.537083E-02/ DATA (BRNDM(i),i=1,64) / & 0.789571E-01, 0.920256E-01, 0.696960E-01, 0.245544E+00, & 0.188503E+00, 0.266127E+00, 0.271371E+00, 0.330917E+00, & 0.190424E+00, 0.224498E+00, 0.282517E+00, 0.130675E+00, & 0.212579E+00, 0.227298E+00, 0.138585E+00, 0.187106E+00, & 0.194527E+00, 0.177034E+00, 0.115902E+00, 0.118499E+00, & 0.142848E+00, 0.216869E+00, 0.149848E+00, 0.971585E-01, & 0.151532E+00, 0.865628E-01, 0.764246E-01, 0.100035E+00, & 0.171133E+00, 0.134737E+00, 0.105173E+00, 0.860832E-01, & 0.148921E+00, 0.869234E-01, 0.106018E+00, 0.184865E+00, & 0.767454E-01, 0.108981E+00, 0.123094E+00, 0.177287E+00, & 0.848146E-01, 0.119356E+00, 0.133829E+00, 0.954505E-01, & 0.155405E+00, 0.164167E+00, 0.161390E+00, 0.113287E+00, & 0.714720E-01, 0.741598E-01, 0.719590E-01, 0.140616E+00, & 0.355356E-01, 0.832779E-01, 0.128680E+00, 0.983013E-01, & 0.629660E-01, 0.643346E-01, 0.717082E-01, 0.629730E-01, & 0.875182E-01, 0.857907E-01, 0.358808E+00, 0.178840E+00/ DATA (BRNDM(i),i=65,128) / & 0.254265E+00, 0.297901E+00, 0.153916E+00, 0.537774E+00, & 0.267906E+00, 0.104254E+00, 0.400723E+00, 0.389670E+00, & 0.263701E+00, 0.338116E+00, 0.351528E+00, 0.267764E+00, & 0.186419E+00, 0.238237E+00, 0.210408E+00, 0.176869E+00, & 0.114715E+00, 0.173299E+00, 0.967770E-01, 0.172565E+00, & 0.162085E+00, 0.157782E+00, 0.886832E-01, 0.242999E+00, & 0.760298E-01, 0.164248E+00, 0.221428E+00, 0.166799E+00, & 0.312514E+00, 0.380600E+00, 0.353828E+00, 0.269500E+00, & 0.254759E+00, 0.285408E+00, 0.159764E+00, 0.721058E-01, & 0.170528E+00, 0.231595E+00, 0.307184E+00, 0.564136E-01, & 0.159884E+00, 0.147907E+00, 0.185666E+00, 0.183567E+00, & 0.182482E+00, 0.230650E+00, 0.175348E+00, 0.195978E+00, & 0.255323E+00, 0.198517E+00, 0.195500E+00, 0.208356E+00, & 0.309603E+00, 0.112011E+00, 0.102570E+00, 0.128276E+00, & 0.168100E+00, 0.177836E+00, 0.105533E+00, 0.903330E-01, & 0.126036E+00, 0.101430E+00, 0.124546E+00, 0.221406E+00/ DATA (BRNDM(i),i=129,163) / & 0.137509E+00, 0.911365E-01, 0.724508E-01, 0.795788E-01, & 0.137411E+00, 0.549175E-01, 0.787714E-01, 0.165544E+00, & 0.136484E+00, 0.146729E+00, 0.820496E-01, 0.846211E-01, & 0.785821E-01, 0.122527E+00, 0.125359E+00, 0.101589E+00, & 0.155756E+00, 0.189239E+00, 0.999086E-01, 0.480993E+00, & 0.100233E+00, 0.153754E+00, 0.130780E+00, 0.136136E+00, & 0.159353E+00, 0.156634E+00, 0.272265E+00, 0.186874E+00, & 0.192090E+00, 0.135397E+00, 0.131497E+00, 0.127463E+00, & 0.227233E+00, 0.190562E+00, 0.214005E+00/ DATA (BANDLO(i),i=1,64) / & 0.000000E+00, 0.100000E+02, 0.200000E+02, 0.300000E+02, & 0.400000E+02, 0.500000E+02, 0.600000E+02, 0.700000E+02, & 0.800000E+02, 0.900000E+02, 0.100000E+03, 0.110000E+03, & 0.120000E+03, 0.130000E+03, 0.140000E+03, 0.150000E+03, & 0.160000E+03, 0.170000E+03, 0.180000E+03, 0.190000E+03, & 0.200000E+03, 0.210000E+03, 0.220000E+03, 0.230000E+03, & 0.240000E+03, 0.250000E+03, 0.260000E+03, 0.270000E+03, & 0.280000E+03, 0.290000E+03, 0.300000E+03, 0.310000E+03, & 0.320000E+03, 0.330000E+03, 0.340000E+03, 0.350000E+03, & 0.360000E+03, 0.370000E+03, 0.380000E+03, 0.390000E+03, & 0.400000E+03, 0.410000E+03, 0.420000E+03, 0.430000E+03, & 0.440000E+03, 0.450000E+03, 0.460000E+03, 0.470000E+03, & 0.480000E+03, 0.490000E+03, 0.500000E+03, 0.510000E+03, & 0.520000E+03, 0.530000E+03, 0.540000E+03, 0.550000E+03, & 0.560000E+03, 0.670000E+03, 0.800000E+03, 0.900000E+03, & 0.990000E+03, 0.107000E+04, 0.120000E+04, 0.121000E+04/ DATA (BANDLO(i),i=65,128) / & 0.122000E+04, 0.123000E+04, 0.124000E+04, 0.125000E+04, & 0.126000E+04, 0.127000E+04, 0.128000E+04, 0.129000E+04, & 0.130000E+04, 0.131000E+04, 0.132000E+04, 0.133000E+04, & 0.134000E+04, 0.135000E+04, 0.136000E+04, 0.137000E+04, & 0.138000E+04, 0.139000E+04, 0.140000E+04, 0.141000E+04, & 0.142000E+04, 0.143000E+04, 0.144000E+04, 0.145000E+04, & 0.146000E+04, 0.147000E+04, 0.148000E+04, 0.149000E+04, & 0.150000E+04, 0.151000E+04, 0.152000E+04, 0.153000E+04, & 0.154000E+04, 0.155000E+04, 0.156000E+04, 0.157000E+04, & 0.158000E+04, 0.159000E+04, 0.160000E+04, 0.161000E+04, & 0.162000E+04, 0.163000E+04, 0.164000E+04, 0.165000E+04, & 0.166000E+04, 0.167000E+04, 0.168000E+04, 0.169000E+04, & 0.170000E+04, 0.171000E+04, 0.172000E+04, 0.173000E+04, & 0.174000E+04, 0.175000E+04, 0.176000E+04, 0.177000E+04, & 0.178000E+04, 0.179000E+04, 0.180000E+04, 0.181000E+04, & 0.182000E+04, 0.183000E+04, 0.184000E+04, 0.185000E+04/ DATA (BANDLO(i),i=129,163) / & 0.186000E+04, 0.187000E+04, 0.188000E+04, 0.189000E+04, & 0.190000E+04, 0.191000E+04, 0.192000E+04, 0.193000E+04, & 0.194000E+04, 0.195000E+04, 0.196000E+04, 0.197000E+04, & 0.198000E+04, 0.199000E+04, 0.200000E+04, 0.201000E+04, & 0.202000E+04, 0.203000E+04, 0.204000E+04, 0.205000E+04, & 0.206000E+04, 0.207000E+04, 0.208000E+04, 0.209000E+04, & 0.210000E+04, 0.211000E+04, 0.212000E+04, 0.213000E+04, & 0.214000E+04, 0.215000E+04, 0.216000E+04, 0.217000E+04, & 0.218000E+04, 0.219000E+04, 0.227000E+04/ DATA (BANDHI(i),i=1,64) / & 0.100000E+02, 0.200000E+02, 0.300000E+02, 0.400000E+02, & 0.500000E+02, 0.600000E+02, 0.700000E+02, 0.800000E+02, & 0.900000E+02, 0.100000E+03, 0.110000E+03, 0.120000E+03, & 0.130000E+03, 0.140000E+03, 0.150000E+03, 0.160000E+03, & 0.170000E+03, 0.180000E+03, 0.190000E+03, 0.200000E+03, & 0.210000E+03, 0.220000E+03, 0.230000E+03, 0.240000E+03, & 0.250000E+03, 0.260000E+03, 0.270000E+03, 0.280000E+03, & 0.290000E+03, 0.300000E+03, 0.310000E+03, 0.320000E+03, & 0.330000E+03, 0.340000E+03, 0.350000E+03, 0.360000E+03, & 0.370000E+03, 0.380000E+03, 0.390000E+03, 0.400000E+03, & 0.410000E+03, 0.420000E+03, 0.430000E+03, 0.440000E+03, & 0.450000E+03, 0.460000E+03, 0.470000E+03, 0.480000E+03, & 0.490000E+03, 0.500000E+03, 0.510000E+03, 0.520000E+03, & 0.530000E+03, 0.540000E+03, 0.550000E+03, 0.560000E+03, & 0.670000E+03, 0.800000E+03, 0.900000E+03, 0.990000E+03, & 0.107000E+04, 0.120000E+04, 0.121000E+04, 0.122000E+04/ DATA (BANDHI(i),i=65,128) / & 0.123000E+04, 0.124000E+04, 0.125000E+04, 0.126000E+04, & 0.127000E+04, 0.128000E+04, 0.129000E+04, 0.130000E+04, & 0.131000E+04, 0.132000E+04, 0.133000E+04, 0.134000E+04, & 0.135000E+04, 0.136000E+04, 0.137000E+04, 0.138000E+04, & 0.139000E+04, 0.140000E+04, 0.141000E+04, 0.142000E+04, & 0.143000E+04, 0.144000E+04, 0.145000E+04, 0.146000E+04, & 0.147000E+04, 0.148000E+04, 0.149000E+04, 0.150000E+04, & 0.151000E+04, 0.152000E+04, 0.153000E+04, 0.154000E+04, & 0.155000E+04, 0.156000E+04, 0.157000E+04, 0.158000E+04, & 0.159000E+04, 0.160000E+04, 0.161000E+04, 0.162000E+04, & 0.163000E+04, 0.164000E+04, 0.165000E+04, 0.166000E+04, & 0.167000E+04, 0.168000E+04, 0.169000E+04, 0.170000E+04, & 0.171000E+04, 0.172000E+04, 0.173000E+04, 0.174000E+04, & 0.175000E+04, 0.176000E+04, 0.177000E+04, 0.178000E+04, & 0.179000E+04, 0.180000E+04, 0.181000E+04, 0.182000E+04, & 0.183000E+04, 0.184000E+04, 0.185000E+04, 0.186000E+04/ DATA (BANDHI(i),i=129,163) / & 0.187000E+04, 0.188000E+04, 0.189000E+04, 0.190000E+04, & 0.191000E+04, 0.192000E+04, 0.193000E+04, 0.194000E+04, & 0.195000E+04, 0.196000E+04, 0.197000E+04, 0.198000E+04, & 0.199000E+04, 0.200000E+04, 0.201000E+04, 0.202000E+04, & 0.203000E+04, 0.204000E+04, 0.205000E+04, 0.206000E+04, & 0.207000E+04, 0.208000E+04, 0.209000E+04, 0.210000E+04, & 0.211000E+04, 0.212000E+04, 0.213000E+04, 0.214000E+04, & 0.215000E+04, 0.216000E+04, 0.217000E+04, 0.218000E+04, & 0.219000E+04, 0.220000E+04, 0.238000E+04/ DATA (AP(i),i=1,64) / & -0.675950E-02, -0.909459E-02, -0.800214E-02, -0.658673E-02, & -0.245580E-02, -0.710464E-02, -0.205565E-02, -0.446529E-02, & -0.440265E-02, -0.593625E-02, -0.201913E-02, -0.349169E-02, & -0.209324E-02, -0.127980E-02, -0.388007E-02, -0.140542E-02, & 0.518346E-02, -0.159375E-02, 0.250508E-02, 0.132182E-01, & -0.903779E-03, 0.110959E-01, 0.924528E-03, 0.207428E-01, & 0.364166E-02, 0.365229E-02, 0.884367E-02, 0.617260E-02, & 0.701340E-02, 0.184265E-01, 0.992822E-02, 0.908582E-02, & 0.106581E-01, 0.276268E-02, 0.158414E-01, 0.145747E-01, & 0.453080E-02, 0.214767E-01, 0.553895E-02, 0.195031E-01, & 0.237016E-01, 0.112371E-01, 0.275977E-01, 0.188833E-01, & 0.131079E-01, 0.130019E-01, 0.385122E-01, 0.111768E-01, & 0.622620E-02, 0.194397E-01, 0.134360E-01, 0.207829E-01, & 0.147960E-01, 0.744479E-02, 0.107564E-01, 0.181562E-01, & 0.170062E-01, 0.233303E-01, 0.256735E-01, 0.274745E-01, & 0.279259E-01, 0.197002E-01, 0.140268E-01, 0.185933E-01/ DATA (AP(i),i=65,128) / & 0.169525E-01, 0.214410E-01, 0.136577E-01, 0.169510E-01, & 0.173025E-01, 0.958346E-02, 0.255024E-01, 0.308943E-01, & 0.196031E-01, 0.183608E-01, 0.149419E-01, 0.206358E-01, & 0.140654E-01, 0.172797E-01, 0.145470E-01, 0.982987E-02, & 0.116695E-01, 0.811333E-02, 0.965823E-02, 0.649977E-02, & 0.462192E-02, 0.545929E-02, 0.680407E-02, 0.291235E-02, & -0.974773E-03, 0.341591E-02, 0.376198E-02, 0.770610E-03, & -0.940864E-04, 0.514532E-02, 0.232371E-02, -0.177741E-02, & -0.374892E-03, -0.370485E-03, -0.221435E-02, -0.490000E-02, & 0.588664E-02, 0.931411E-03, -0.456043E-03, -0.545576E-02, & -0.421136E-02, -0.353742E-02, -0.174276E-02, -0.361246E-02, & -0.337822E-02, -0.867030E-03, -0.118001E-02, -0.222405E-02, & -0.725144E-03, 0.118483E-02, 0.995087E-02, 0.273812E-03, & 0.417298E-02, 0.764294E-02, 0.631568E-02, -0.213528E-02, & 0.746130E-02, 0.110337E-02, 0.153157E-01, 0.504532E-02, & 0.406047E-02, 0.192895E-02, 0.202058E-02, 0.126420E-01/ DATA (AP(i),i=129,163) / & 0.310028E-02, 0.214779E-01, 0.560165E-02, 0.661070E-02, & 0.694966E-02, 0.539194E-02, 0.103745E-01, 0.180150E-01, & 0.747133E-02, 0.114927E-01, 0.115213E-01, 0.160709E-02, & 0.154278E-01, 0.112067E-01, 0.148690E-01, 0.154442E-01, & 0.123977E-01, 0.237539E-01, 0.162820E-01, 0.269484E-01, & 0.178081E-01, 0.143221E-01, 0.262468E-01, 0.217065E-01, & 0.107083E-01, 0.281220E-01, 0.115565E-01, 0.231244E-01, & 0.225197E-01, 0.178624E-01, 0.327708E-01, 0.116657E-01, & 0.277452E-01, 0.301647E-01, 0.349782E-01/ DATA (BP(i),i=1,64) / & 0.717848E-05, 0.169280E-04, 0.126710E-04, 0.758397E-05, & -0.533900E-05, 0.143490E-04, -0.595854E-05, 0.296465E-05, & 0.323446E-05, 0.115359E-04, -0.692861E-05, 0.131477E-04, & -0.624945E-05, -0.756955E-06, 0.107458E-05, -0.159796E-05, & -0.290529E-04, -0.170918E-05, -0.193934E-04, -0.707209E-04, & -0.148154E-04, -0.383162E-04, -0.186050E-04, -0.951796E-04, & -0.210944E-04, -0.330590E-04, -0.373087E-04, -0.408972E-04, & -0.396759E-04, -0.827756E-04, -0.573773E-04, -0.325384E-04, & -0.449411E-04, -0.271450E-04, -0.752791E-04, -0.549699E-04, & -0.225655E-04, -0.102034E-03, -0.740322E-05, -0.668846E-04, & -0.106063E-03, -0.304840E-04, -0.796023E-04, 0.504880E-04, & 0.486384E-04, -0.531946E-04, -0.147771E-03, -0.406785E-04, & 0.615750E-05, -0.486264E-04, -0.419335E-04, -0.819467E-04, & -0.709498E-04, 0.326984E-05, -0.369743E-04, -0.526848E-04, & -0.550050E-04, -0.684057E-04, -0.447093E-04, -0.778390E-04, & -0.982953E-04, -0.772497E-04, -0.119430E-05, -0.655187E-04/ DATA (BP(i),i=65,128) / & -0.339078E-04, 0.716657E-04, -0.335893E-04, 0.220239E-04, & -0.491012E-04, -0.393325E-04, -0.626461E-04, -0.795479E-04, & -0.599181E-04, -0.578153E-04, -0.597559E-05, -0.866750E-04, & -0.486783E-04, -0.580912E-04, -0.647368E-04, -0.350643E-04, & -0.566635E-04, -0.385738E-04, -0.463782E-04, -0.321485E-04, & -0.177300E-04, -0.250201E-04, -0.365492E-04, -0.165218E-04, & -0.649177E-05, -0.218458E-04, -0.984604E-05, -0.120034E-04, & -0.110119E-06, -0.164405E-04, -0.141396E-04, 0.315347E-05, & -0.141544E-05, -0.297320E-05, -0.216248E-05, 0.839264E-05, & -0.178197E-04, -0.106225E-04, -0.468195E-05, 0.997043E-05, & 0.679709E-05, 0.324610E-05, -0.367325E-05, 0.671058E-05, & 0.509293E-05, -0.437392E-05, -0.787922E-06, -0.271503E-06, & -0.437940E-05, -0.128205E-04, -0.417830E-04, -0.561134E-05, & -0.209940E-04, -0.414366E-04, -0.289765E-04, 0.680406E-06, & -0.558644E-05, -0.530395E-05, -0.622242E-04, -0.159979E-05, & -0.140286E-04, -0.128463E-04, -0.929499E-05, -0.327886E-04/ DATA (BP(i),i=129,163) / & -0.189353E-04, -0.737589E-04, -0.323471E-04, -0.272502E-04, & -0.321731E-04, -0.326958E-04, -0.509157E-04, -0.681890E-04, & -0.362182E-04, -0.354405E-04, -0.578392E-04, 0.238627E-05, & -0.709028E-04, -0.518717E-04, -0.491859E-04, -0.718017E-04, & -0.418978E-05, -0.940819E-04, -0.630375E-04, -0.478469E-04, & -0.751896E-04, -0.267113E-04, -0.109019E-03, -0.890983E-04, & -0.177301E-04, -0.120216E-03, 0.220464E-04, -0.734277E-04, & -0.868068E-04, -0.652319E-04, -0.136982E-03, -0.279933E-06, & -0.791824E-04, -0.111781E-03, -0.748263E-04/ DATA (ATP(i),i=1,64) / & -0.722782E-02, -0.901531E-02, -0.821263E-02, -0.808024E-02, & -0.320169E-02, -0.661305E-02, -0.287272E-02, -0.486143E-02, & -0.242857E-02, -0.530288E-02, -0.146813E-02, -0.566474E-03, & -0.102192E-02, 0.300643E-03, -0.331655E-02, 0.648220E-03, & 0.552446E-02, -0.933046E-03, 0.205703E-02, 0.130638E-01, & -0.229828E-02, 0.715648E-02, 0.444446E-03, 0.193500E-01, & 0.364119E-02, 0.252713E-02, 0.102420E-01, 0.494224E-02, & 0.584934E-02, 0.146255E-01, 0.921986E-02, 0.768012E-02, & 0.916105E-02, 0.276223E-02, 0.125245E-01, 0.131146E-01, & 0.793016E-02, 0.201536E-01, 0.658631E-02, 0.171711E-01, & 0.228470E-01, 0.131306E-01, 0.226658E-01, 0.176086E-01, & 0.149987E-01, 0.143060E-01, 0.313189E-01, 0.117070E-01, & 0.133522E-01, 0.244259E-01, 0.148393E-01, 0.223982E-01, & 0.151792E-01, 0.180474E-01, 0.106299E-01, 0.191016E-01, & 0.171776E-01, 0.229724E-01, 0.275530E-01, 0.302731E-01, & 0.281662E-01, 0.199525E-01, 0.192588E-01, 0.173220E-01/ DATA (ATP(i),i=65,128) / & 0.195220E-01, 0.169371E-01, 0.193212E-01, 0.145558E-01, & 0.189654E-01, 0.122030E-01, 0.186206E-01, 0.228842E-01, & 0.139343E-01, 0.164006E-01, 0.137276E-01, 0.154005E-01, & 0.114575E-01, 0.129956E-01, 0.115305E-01, 0.929260E-02, & 0.106359E-01, 0.771623E-02, 0.106075E-01, 0.597630E-02, & 0.493960E-02, 0.532554E-02, 0.646175E-02, 0.302693E-02, & 0.150899E-02, 0.310333E-02, 0.533734E-02, 0.239094E-03, & 0.356782E-02, 0.707574E-02, 0.215758E-02, -0.527589E-03, & 0.643893E-03, -0.101916E-02, -0.383336E-02, -0.445966E-02, & 0.880190E-02, 0.245662E-02, -0.560923E-03, -0.582201E-02, & -0.323233E-02, -0.454197E-02, -0.240905E-02, -0.343160E-02, & -0.335156E-02, -0.623846E-03, 0.393633E-03, -0.271593E-02, & -0.675874E-03, 0.920642E-03, 0.102168E-01, -0.250663E-03, & 0.437126E-02, 0.767434E-02, 0.569931E-02, -0.929326E-03, & 0.659414E-02, 0.280687E-02, 0.127614E-01, 0.780789E-02, & 0.374807E-02, 0.274288E-02, 0.534940E-02, 0.104349E-01/ DATA (ATP(i),i=129,163) / & 0.294379E-02, 0.177846E-01, 0.523249E-02, 0.125339E-01, & 0.548538E-02, 0.577403E-02, 0.101532E-01, 0.170375E-01, & 0.758396E-02, 0.113402E-01, 0.106960E-01, 0.107782E-01, & 0.136148E-01, 0.992064E-02, 0.167276E-01, 0.149603E-01, & 0.136259E-01, 0.234521E-01, 0.166806E-01, 0.298505E-01, & 0.167592E-01, 0.186679E-01, 0.233062E-01, 0.228467E-01, & 0.128947E-01, 0.293979E-01, 0.219815E-01, 0.220663E-01, & 0.272710E-01, 0.237139E-01, 0.331743E-01, 0.208799E-01, & 0.281472E-01, 0.318440E-01, 0.370962E-01/ DATA (BTP(i),i=1,64) / & 0.149748E-04, 0.188007E-04, 0.196530E-04, 0.124747E-04, & -0.215751E-07, 0.128357E-04, -0.265798E-05, 0.606262E-05, & 0.287668E-05, 0.974612E-05, -0.833451E-05, 0.584410E-05, & -0.452879E-05, -0.782537E-05, 0.786165E-05, -0.768351E-05, & -0.196168E-04, 0.177297E-06, -0.129258E-04, -0.642798E-04, & -0.986297E-05, -0.257145E-04, -0.141996E-04, -0.865089E-04, & -0.141691E-04, -0.272578E-04, -0.295198E-04, -0.308878E-04, & -0.313193E-04, -0.669272E-04, -0.475777E-04, -0.221332E-04, & -0.419930E-04, -0.102519E-04, -0.590184E-04, -0.574771E-04, & -0.240809E-04, -0.913994E-04, -0.908886E-05, -0.721074E-04, & -0.902837E-04, -0.447582E-04, -0.664544E-04, -0.143150E-04, & -0.511866E-05, -0.559352E-04, -0.104734E-03, -0.305206E-04, & 0.103303E-04, -0.613019E-04, -0.320040E-04, -0.738909E-04, & -0.388263E-04, 0.306515E-04, -0.352214E-04, -0.253940E-04, & -0.521369E-04, -0.746260E-04, -0.744124E-04, -0.881905E-04, & -0.933645E-04, -0.664045E-04, -0.570712E-05, -0.566312E-04/ DATA (BTP(i),i=65,128) / & -0.364967E-04, 0.393501E-06, -0.234050E-04, -0.141317E-04, & -0.525480E-04, -0.172241E-04, -0.410843E-04, -0.358348E-04, & -0.256168E-04, -0.509482E-04, -0.180570E-04, -0.555356E-04, & -0.271464E-04, -0.274040E-04, -0.480889E-04, -0.275751E-04, & -0.415681E-04, -0.383770E-04, -0.280139E-04, -0.287919E-04, & -0.125865E-04, -0.265467E-04, -0.172765E-04, -0.164611E-04, & 0.189183E-04, -0.171219E-04, -0.132766E-04, -0.344611E-05, & -0.442832E-05, -0.185779E-04, -0.139755E-04, 0.168083E-05, & -0.395287E-05, -0.297871E-05, 0.434383E-05, 0.131741E-04, & -0.192637E-04, -0.549551E-05, 0.122553E-05, 0.204627E-04, & 0.154027E-04, 0.953462E-05, 0.131125E-05, 0.732839E-05, & 0.755405E-05, -0.305552E-05, -0.434858E-05, 0.308409E-05, & -0.164787E-05, -0.818533E-05, -0.355041E-04, -0.504696E-05, & -0.229022E-04, -0.356891E-04, -0.230346E-04, 0.518835E-05, & -0.160187E-04, -0.104617E-04, -0.464754E-04, -0.115807E-04, & -0.130230E-04, -0.603491E-05, -0.125324E-04, -0.165516E-04/ DATA (BTP(i),i=129,163) / & -0.991679E-05, -0.529432E-04, -0.200199E-04, -0.181977E-04, & -0.220940E-04, -0.204483E-04, -0.432584E-04, -0.449109E-04, & -0.247305E-04, -0.174253E-04, -0.484446E-04, 0.354150E-04, & -0.425581E-04, -0.406562E-04, -0.505495E-04, -0.651856E-04, & -0.153953E-04, -0.894294E-04, -0.616551E-04, -0.846504E-04, & -0.699414E-04, -0.376203E-04, -0.940985E-04, -0.753050E-04, & -0.183710E-04, -0.123907E-03, -0.279347E-04, -0.736381E-04, & -0.103588E-03, -0.754117E-04, -0.140991E-03, -0.366687E-04, & -0.927785E-04, -0.125321E-03, -0.115290E-03/ DATA (BETAD(i),i=1,64) / & 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.000000E+00, & 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.000000E+00, & 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.000000E+00, & 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.000000E+00, & 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.000000E+00, & 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.000000E+00, & 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.000000E+00, & 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.000000E+00, & 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.000000E+00, & 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.000000E+00, & 0.234879E+03, 0.217419E+03, 0.201281E+03, 0.186364E+03, & 0.172576E+03, 0.159831E+03, 0.148051E+03, 0.137163E+03, & 0.127099E+03, 0.117796E+03, 0.109197E+03, 0.101249E+03, & 0.939031E+02, 0.871127E+02, 0.808363E+02, 0.750349E+02, & 0.497489E+02, 0.221212E+02, 0.113124E+02, 0.754174E+01, & 0.589554E+01, 0.495227E+01, 0.000000E+00, 0.000000E+00/ DATA (BETAD(i),i=65,128) / & 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.000000E+00, & 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.000000E+00, & 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.000000E+00, & 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.000000E+00, & 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.000000E+00, & 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.000000E+00, & 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.000000E+00, & 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.000000E+00, & 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.000000E+00, & 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.000000E+00, & 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.000000E+00, & 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.000000E+00, & 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.000000E+00, & 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.000000E+00, & 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.000000E+00, & 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.000000E+00/ DATA (BETAD(i),i=129,163) / & 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.000000E+00, & 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.000000E+00, & 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.000000E+00, & 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.000000E+00, & 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.000000E+00, & 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.000000E+00, & 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.000000E+00, & 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.000000E+00, & 0.000000E+00, 0.000000E+00, 0.000000E+00/ !----------------------------------------------------------------------- ! ! Random band parameters for the longwave calculations using ! comboned wide frequency bands between 160 and 1200 cm-1, ! as well as the 2270-2380 band for source calculations. ! ! BANDS 1-8: COMBINED WIDE FREQUENCY BANDS FOR 160-560 CM-1 ! BANDS 9-14: FREQUENCY BANDS,AS IN BANDTA (NARROW BANDS) ! FOR 560-1200 CM-1 ! BAND 15: FREQUENCY BAND 2270-2380 CM-1,USED FOR SOURCE ! CALCULATION ONLY ! ! Thus NBLY presently equals 15 ! ! Bands are arranged in order of increasing wavenumber. ! ! ACOMB = RANDOM "A" PARAMETER FOR (NBLY) BANDS ! BCOMB = RANDOM "B" PARAMETER FOR (NBLY) BANDS ! BETACM = CONTINUUM COEFFICIENTS FOR (NBLY) BANDS ! APCM,BPCM = CAPPHI COEFFICIENTS FOR (NBLY) BANDS ! ATPCM,BTPCM = CAPPSI COEFFICIENTS FOR (NBLY) BANDS ! BDLOCM = LOWEST FREQUENCY IN EACH OF (NBLY) FREQ. BANDS ! BDHICM = HIGHEST FREQUENCY IN EACH OF (NBLY) FREQ. BANDS ! AO3CM = RANDOM "A" PARAMETER FOR OZONE IN (3) OZONE ! BANDS. ! BO3CM = RANDOM "B" PARAMETER FOR OZONE IN (3) OZONE ! BANDS ! AB15CM = THE PRODUCT ARNDM*BRNDM FOR THE TWO BANDS ! REPRESENTING THE 15 UM BAND COMPLEX OF CO2 ! BETINC = CONT.COEFFICIENT FOR A SPECIFIED WIDE ! FREQ.BAND (800-990 AND 1070-1200 CM-1). ! IBAND = INDEX NO OF THE 40 WIDE BANDS USED IN ! COMBINED WIDE BAND CALCULATIONS. IN OTHER ! WORDS,INDEX TELLING WHICH OF THE 40 WIDE ! BANDS BETWEEN 160-560 CM-1 ARE INCLUDED IN ! EACH OF THE FIRST 8 COMBINED WIDE BANDS ! ! Data for ACOMB,BCOMB,APCM,BPCM,ATPCM,BTPCM,AO3CM,BO3CM are ! obtained by using the AFGL 1982 catalog. Continuum coefficients ! are from Roberts (1976). IBAND index values are obtained by ! experimentation. This data was formerly in COMMON /BDCOMB/. INTEGER IBAND(40) REAL ACOMB(NBLY),BCOMB(NBLY), & BETACM(NBLY),APCM(NBLY),BPCM(NBLY),ATPCM(NBLY), & BTPCM(NBLY),BDLOCM(NBLY),BDHICM(NBLY),BETINC, & AO3CM(3),BO3CM(3),AB15CM(2) DATA ACOMB / & 0.152070E+05, 0.332194E+04, 0.527177E+03, 0.163124E+03, & 0.268808E+03, 0.534591E+02, 0.268071E+02, 0.123133E+02, & 0.600199E+01, 0.640803E+00, 0.501549E-01, 0.167961E-01, & 0.178110E-01, 0.170166E+00, 0.537083E-02/ DATA BCOMB / & 0.152538E+00, 0.118677E+00, 0.103660E+00, 0.100119E+00, & 0.127518E+00, 0.118409E+00, 0.904061E-01, 0.642011E-01, & 0.629660E-01, 0.643346E-01, 0.717082E-01, 0.629730E-01, & 0.875182E-01, 0.857907E-01, 0.214005E+00/ DATA APCM / & -0.671879E-03, 0.654345E-02, 0.143657E-01, 0.923593E-02, & 0.117022E-01, 0.159596E-01, 0.181600E-01, 0.145013E-01, & 0.170062E-01, 0.233303E-01, 0.256735E-01, 0.274745E-01, & 0.279259E-01, 0.197002E-01, 0.349782E-01/ DATA BPCM / & -0.113520E-04, -0.323965E-04, -0.448417E-04, -0.230779E-04, & -0.361981E-04, -0.145117E-04, 0.198349E-04, -0.486529E-04, & -0.550050E-04, -0.684057E-04, -0.447093E-04, -0.778390E-04, & -0.982953E-04, -0.772497E-04, -0.748263E-04/ DATA ATPCM / & -0.106346E-02, 0.641531E-02, 0.137362E-01, 0.922513E-02, & 0.136162E-01, 0.169791E-01, 0.206959E-01, 0.166223E-01, & 0.171776E-01, 0.229724E-01, 0.275530E-01, 0.302731E-01, & 0.281662E-01, 0.199525E-01, 0.370962E-01/ DATA BTPCM / & -0.735731E-05, -0.294149E-04, -0.505592E-04, -0.280894E-04, & -0.492972E-04, -0.341508E-04, -0.362947E-04, -0.250487E-04, & -0.521369E-04, -0.746260E-04, -0.744124E-04, -0.881905E-04, & -0.933645E-04, -0.664045E-04, -0.115290E-03/ DATA BETACM / & 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.000000E+00, & 0.188625E+03, 0.144293E+03, 0.174098E+03, 0.909366E+02, & 0.497489E+02, 0.221212E+02, 0.113124E+02, 0.754174E+01, & 0.589554E+01, 0.495227E+01, 0.000000E+00/ DATA IBAND / & 2, 1, 2, 2, 1, 2, 1, 3, 2, 2, & 3, 2, 2, 4, 2, 4, 2, 3, 3, 2, & 4, 3, 4, 3, 7, 5, 6, 7, 6, 5, & 7, 6, 7, 8, 6, 6, 8, 8, 8, 8/ !----------------------------------------------------------------------- ! ! Random band parameters for specific wide bands. At present, ! the information consists of: 1) random model parameters for ! the 15 um band, 560-800 cm-1; 2) the continuum coefficient for ! the 800-990, 1070-1200 cm-1 band. ! ! specifically: ! AWIDE = RANDOM "A" PARAMETER FOR BAND ! BWIDE = RANDOM "B" PARAMETER FOR BAND ! BETAWD = CONTINUUM COEFFICIENTS FOR BAND ! APWD,BPWD = CAPPHI COEFFICIENTS FOR BAND ! ATPWD,BTPWD = CAPPSI COEFFICIENTS FOR BAND ! BDLOWD = LOWEST FREQUENCY IN EACH FREQ BAND ! BDHIWD = HIGHEST FREQUENCY IN EACH FREQ BAND ! AB15WD = THE PRODUCT ARNDM*BRNDM FOR THE ONE BAND ! REPRESENTING THE 15 UM BAND COMPLEX OF CO2 ! BETINW = CONT.COEFFICIENT FOR A SPECIFIED WIDE ! FREQ.BAND (800-990 AND 1070-1200 CM-1). ! SKO2D = 1./BETINW, USED IN SPA88 FOR CONT. COEFFS ! SKC1R = BETAWD/BETINW, USED FOR CONT. COEFF. FOR ! 15 UM BAND IN FST88 ! SKO3R = RATIO OF CONT. COEFF. FOR 9.9 UM BAND TO ! BETINW, USED FOR 9.6 UM CONT COEFF IN FST88 ! ! Data for AWIDE,BWIDE,APWD,BPWD,ATPWD,BTPWD,AO3WD,BO3WD are ! obtained by using the AFGL 1982 catalog. Continuum coefficients ! are from Roberts (1976). This data was formerly in ! COMMON /BDWIDE/. REAL AWIDE,BWIDE,BETAWD, & APWD,BPWD,ATPWD,BTPWD, & BDLOWD,BDHIWD,BETINW, & AB15WD,SKO2D,SKC1R,SKO3R DATA AWIDE / 0.309801E+01/ DATA BWIDE / 0.495357E-01/ DATA APWD / 0.177115E-01/ DATA BPWD /-0.545226E-04/ DATA ATPWD / 0.187967E-01/ DATA BTPWD /-0.567449E-04/ DATA BETAWD / 0.347839E+02/ DATA BETINW / 0.766811E+01/ DATA BDLOWD / 0.560000E+03/ DATA BDHIWD / 0.800000E+03/ !----------------------------------------------------------------------- ! ! CLDFAC = CLOUD TRANSMISSION FUNCTION,ASSUMING RANDOM ! OVERLAP REAL, ALLOCATABLE, DIMENSION(:,:,:) :: CLDFAC !----------------------------------------------------------------------- ! ! Basic quantities computed in SUBROUTINE LWRAD and used in ! the remaining longwave routines (formally COMMON /KDACOM/): ! ! QH2O = H2O MASS MIXING RATIO,MULTIPLIED BY THE ! DIFFUSIVITY FACTOR (DIFFCTR) ! P = PRESSURE AT FLUX LEVELS OF MODEL ! DELP2 = PRESSURE DIFFERENCE BETWEEN FLUX LEVELS ! DELP = INVERSE OF DELP2 ! TTTT = TEMPERATURE ASSIGNED TO MODEL FLUX LEVELS ! VAR1 = H2O OPTICAL PATH IN MODEL LAYERS (BETWEEN ! FLUX LEVELS) ! VAR2 = PRESSURE-WEIGHTED H2O OPTICAL PATH IN MODEL LAYERS ! VAR3 = O3 OPTICAL PATH IN MODEL LAYERS ! VAR4 = PRESSURE-WEIGHTED O3 OPTICAL PATH IN MODEL LAYERS ! CNTVAL = H2O CONTINUUM PATH IN MODEL LAYERS FOR THE ! 800-990 AND 1070-1200 CM-1 COMBINED BAND REAL, ALLOCATABLE, DIMENSION(:,:) :: QH2O,P,DELP2,DELP,TTTT REAL, ALLOCATABLE, DIMENSION(:,:) :: VAR1,VAR2,VAR3,VAR4,CNTVAL !----------------------------------------------------------------------- ! ! Flux quantities computed by the radiation code, used for ! diagnostic purposes (formally COMMON /RDFLUX/): ! ! FLX1E1 = FLUX AT TOP FOR 0-160,1200-2200 CM-1 RANGE ! GXCTS = FLUX AT TOP FOR 160-1200 CM-1 RANGE ! FCTSG = CTS FLUX AT GROUND. USED TO OBTAIN GXCTS ! BY BANDS. REAL, ALLOCATABLE, DIMENSION(:) :: FLX1E1,GXCTS REAL, ALLOCATABLE, DIMENSION(:,:) :: FCTSG !----------------------------------------------------------------------- ! ! Planck function values used for the radiative calculations ! (formally COMMON /SRCCOM/): ! ! SORC = PLANCK FCTN, AT MODEL TEMPERATURES, FOR ALL BANDS ! USED IN CTS CALCULATIONS ! CSOUR1 = PLANCK FCTN FOR 560-670 CM-1 BAND ! CSOUR2 = PLANCK FCTN FOR 670-800 CM-1 BAND ! CSOUR = PLANCK FCTN FOR 560-800 CM-1 BANDS ! OSOUR = PLANCK FCTN FOR 990-1070 CM-1 BAND ! SS1 = PLANCK FCTN FOR 800-990,1070-1200 CM-1 BANDS REAL, ALLOCATABLE, DIMENSION(:,:,:) :: SORC REAL, ALLOCATABLE, DIMENSION(:,:) :: CSOUR1,CSOUR2,OSOUR,CSOUR,SS1 !----------------------------------------------------------------------- ! ! Quantities precomputed in subroutine TABLE for use in ! the longwave radiation module (formally COMMON /TABCOM/): ! ! EM1 = E1 FUNCTION, EVALUATED OVER THE 0-560 AND ! 1200-2200 CM-1 INTERVALS ! EM1WDE = E1 FUNCTION, EVALUATED OVER THE 160-560 CM-1 ! INTERVAL ! TABLE1 = E2 FUNCTION, EVALUATED OVER THE 0-560 AND ! 1200-2200 CM-1 INTERVALS ! TABLE2 = TEMPERATURE DERIVATIVE OF TABLE1 ! TABLE3 = MASS DERIVATIVE OF TABLE1 ! EM3 = E3 FUNCTION, EVALUATED OVER THE 0-560 AND ! 1200-2200 CM-1 INTERVALS ! SOURCE = PLANCK FUNCTION, EVALUATED AT SPECIFIED TEMPS. FOR ! BANDS USED IN CTS CALCULATIONS ! DSRCE = TEMPERATURE DERIVATIVE OF SOURCE ! INDX2 = INDEX VALUES USED IN OBTAINING "LOWER TRIANGLE" ! ELEMENTS OF AVEPHI,ETC.,IN FST88 ! KMAXV = INDEX VALUES USED IN OBTAINING "UPPER TRIANGLE" ! ELEMENTS OF AVEPHI,ETC.,IN FST88 ! KMAXVM = KMAXV(LMAX),USED FOR DO LOOP INDICES INTEGER, ALLOCATABLE, DIMENSION(:) :: INDX1,INDX2,KMAXV INTEGER :: KMAXVM INTEGER, PRIVATE :: IMAXold=0,LMAXold=0 INTEGER, PRIVATE :: IMAX REAL, DIMENSION(28,NBLY) :: SOURCE,DSRCE ! REAL EM1 (28,180),EM1WDE(28,180),TABLE1(28,180), ! COMMON /TABCOM/ EM1 (28,180),EM1WDE(28,180),TABLE1(28,180), ! & TABLE2(28,180),TABLE3(28,180),EM3 (28,180) !----------------------------------------------------------------------- ! ! Transmission functions used for radiative computations, and ! output heating rates and fluxes, except those needed out of ! the radiative module (formally COMMON /TFCOM/): ! ! TO3 = TRANSMISSION FCTN FOR THE 990-1070 CM-1 BAND ! O3(9.6 UM) + H2O CONTINUUM (NO LINES) ! CO21 = TRANSMISSION FCTN FOR THE 560-800 CM-1 BAND ! (AS 1 BAND). INCLUDES CO2 (IN LWRAD) AND ! H2O(L+C) AFTER MULTIPLICATION WITH "OVER" ! IN FST88 ! EMISS = E2 EMISSIVITYY FCTN FOR H2O LINES (0-560,1200-2200 ! CM-1). OBTAINED IN E1E288. ! EMISS2 = TRANSMISSION FCTN FOR H2O CONTINUUM IN THE 800-990 ! AND 1070-1200 CM-1 REGION, TAKEN AS 1 BAND ! AVEPHI = H2O OPTICAL PATHS BET. FLUX PRESSURES: INPUT TO ! EMISSIVITY CALCULATIONS. ! TTEMP = TEMPERATURES USED AS INPUT FOR EMISSIVITY CALCS. ! CTS = APPROX CTS HEATING RATES FOR 160-560 AND 800-990, ! 1070-1200 CM-1 RANGES ! CTSO3 = APPROX CTS HEATING RATES FOR 560-800,990-1070 CM-1 ! RANGES ! EXCTS = EXACT CTS HEATING RATES FOR 160-1200 CM-1 RANGE ! EXCTSN = EXACT CTS HEATING RATES, BY BANDS ! E1FLX = E1 EMISSIVITY FCTN FOR H2O LINES (0-560,1200-CM-1) ! CO2NBL = CO2 TRANS. FCTNS. (NOT PRESSURE-INTEGRATED) FOR ! ADJACENT LEVELS,OVER THE 560-800 CM-1 RANGE. ! CO2SP1 = CO2 TRANS. FCTNS. (NOT PRESSURE-INTEGRATED) BET. ! A FLUX LEVEL AND SPACE, FOR THE 560-670 CM-1 ! RANGE. USED FOR EXACT CTS CALCS. ! CO2SP2 = SAME AS CO2SP1, BUT FOR THE 670-800 CM-1 RANGE. ! CO2SP = SAME AS CO2SP1, BUT FOR THE 560-800 CM-1 BAND. ! USED FOR APPROX CTS CALCS. ! TO3SPC = O3 OPTICAL DEPTHS BET. A LEVEL AND SPACE. USED FOR ! EXACT CTS CALCS. ! TOTVO2 = H2O CONTINUUM OPTICAL PATHS BET. SPACE AND A ! LEVEL, USING THE CNT. COEFFICIENT FOR THE ! 1-BAND 800-990,1070-1200 CM-1 BAND. USED FOR ! CTS CALCS. REAL,ALLOCATABLE,DIMENSION(:,:,:) :: TO3,CO21,EMISS,EMISS2,AVEPHI REAL,ALLOCATABLE,DIMENSION(:,:) :: CTS,CTSO3,EXCTS REAL,ALLOCATABLE,DIMENSION(:,:,:) :: EXCTSN REAL,ALLOCATABLE,DIMENSION(:,:) :: E1FLX,CO2NBL,CO2SP1,CO2SP2 REAL,ALLOCATABLE,DIMENSION(:,:) :: CO2SP,TO3SPC,TOTVO2 !------------ VERSION NUMBER ---------------- character(len=128) :: version = '$Id: longwave.F90,v 13.0 2006/03/28 21:09:33 fms Exp $' character(len=128) :: tagname = '$Name: mom4p1_pubrel_dec2009_nnz $' logical :: module_is_initialized = .false. !----------------------------------------------------------------------- public LWRad, Rad_DeAlloc, longwave_init, longwave_end public OSOUR, CSOUR, SS1, FLX1E1, GXCTS, FCTSG, CLDFAC, DELP2, DELP, & TO3, CO21, EMISS, EMISS2, CTS, EXCTS, EXCTSN, E1FLX, CO2SP, & IBAND, BANDLO, BANDHI CONTAINS !####################################################################### !####################################################################### Subroutine longwave_init !------- write version number --------- if ( mpp_pe() == mpp_root_pe() ) then call write_version_number(version, tagname) endif module_is_initialized = .true. !--------------------------------------------------------------------- End Subroutine longwave_init !####################################################################### !####################################################################### Subroutine longwave_end module_is_initialized = .false. !--------------------------------------------------------------------- End Subroutine longwave_end !####################################################################### !####################################################################### SUBROUTINE LWRAD (KTOP,KBTM,NCLDS,EMCLD,PRES,TEMP,RH2O,QO3,CAMT, & RRVCO2, HEATRA,GRNFLX,TOPFLX, LSFC,PSFC) INTEGER, INTENT(IN), DIMENSION(:,:,:) :: KTOP,KBTM INTEGER, INTENT(IN), DIMENSION(:,:) :: NCLDS REAL, INTENT(IN), DIMENSION(:,:,:) :: EMCLD REAL, INTENT(IN), DIMENSION(:,:,:) :: PRES,TEMP,RH2O,QO3,CAMT REAL, INTENT(IN) :: RRVCO2 REAL, INTENT(OUT), DIMENSION(:,:,:) :: HEATRA REAL, INTENT(OUT), DIMENSION(:,:) :: GRNFLX,TOPFLX INTEGER, INTENT(IN),OPTIONAL, DIMENSION(:,:) :: LSFC REAL , INTENT(IN),OPTIONAL, DIMENSION(:,:) :: PSFC INTEGER J IMAX=SIZE(PRES,1) CALL RAD_ALLOC (IMAX) DO j=1,SIZE(PRES,2) CALL CLO88 (KTOP(:,j,:),KBTM(:,j,:),NCLDS(:,j), & CAMT(:,j,:),EMCLD(:,j,:)) IF (PRESENT(LSFC) .and. PRESENT(PSFC)) THEN CALL LWRAD2 (KTOP(:,j,:),KBTM(:,j,:),NCLDS(:,j),EMCLD(:,j,:), & PRES(:,j,:),TEMP(:,j,:),RH2O(:,j,:),QO3(:,j,:), & CAMT(:,j,:), RRVCO2, & HEATRA(:,j,:),GRNFLX(:,j),TOPFLX(:,j), & LSFC(:,j),PSFC(:,j)) ELSE CALL LWRAD2 (KTOP(:,j,:),KBTM(:,j,:),NCLDS(:,j),EMCLD(:,j,:), & PRES(:,j,:),TEMP(:,j,:),RH2O(:,j,:),QO3(:,j,:), & CAMT(:,j,:), RRVCO2, & HEATRA(:,j,:),GRNFLX(:,j),TOPFLX(:,j)) ENDIF ENDDO !del CALL RAD_DEALLOC END SUBROUTINE LWRAD !####################################################################### !####################################################################### SUBROUTINE LWRAD2 (KTOP,KBTM,NCLDS,EMCLD,PRES,TEMP,RH2O,QO3,CAMT, & RRVCO2, HEATRA,GRNFLX,TOPFLX, LSFC,PSFC) !----------------------------------------------------------------------- ! LONG WAVE RADIATION CODE ! ------------------------ ! SUBROUTINE LWRAD (formally LWR88) COMPUTES TEMPERATURE-CORRECTED ! CO2 TRANSMISSION FUNCTIONS AND ALSO COMPUTES THE PRESSURE GRID ! AND LAYER OPTICAL PATHS. !----------------------------------------------------------------------- ! INPUT PARAMETERS ! ---------------- ! ! KTOP = INDEX OF (DATA LEVEL) PRESSURE OF CLOUD TOP,USED ! IN THE LONGWAVE PROGRAM ! KBTM = INDEX OF (DATA LEVEL) PRESSURE OF CLOUD BOTTOM, ! USED IN THE LONGWAVE PROGRAM ! NCLDS = NO. CLOUDS AT EACH GRID PT. ! EMCLD = CLOUD EMISSIVITY. SET TO ONE BY DEFAULT, BUT MAY ! BE MODIFIED FOR USE IN LONGWAVE PROGRAM. ! PRES = PRESSURE (CGS UNITS) AT DATA LEVELS OF MODEL ! TEMP = TEMPERATURE (K) AT DATA LEVELS OF MODEL ! RH2O = MASS MIXING RATIO (G/G) OF H2O AT MODEL DATA LVLS. ! QO3 = MASS MIXING RATIO (G/G) OF O3 AT MODEL DATA LVLS. ! CAMT = CLOUD AMOUNTS OF CLOUDS (THEIR LOCATIONS ARE ! SPECIFIED IN THE KTOP/KBTM INDICES) ! RRVCO2 = THE VOLUME MIXING RATIO OF CO2 (SCALAR) ! (used in NLTE, passed by argument list) ! ! PSFC = Surface pressure, dimensioned by IMAX. ! TSFC = Surface temperature, dimensioned by IMAX. ! LSFC = Vertical index of the lowest model level, ! dimensioned by IMAX. ! !----------------------------------------------------------------------- ! OUTPUT PARAMETERS ! ----------------- ! ! HEATRA = HEATING RATE AT DATA LEVELS (K/DAY) ! GRNFLX = NET LONGWAVE FLUX AT THE GROUND (CGS UNITS) ! TOPFLX = NET LONGWAVE FLUX AT THE TOP (CGS UNITS) ! ! !----------------------------------------------------------------------- ! OTHER INPUTS (FROM COMMON BLOCKS): ! ---------------------------------- ! CO251,CO258,CDT51,CDT58 CO2BD3 ! C2D51,C2D58,CO2M51,CO2M58 CO2BD3 ! CDTM51,CDTM58,C2DM51,C2DM58 CO2BD3 ! STEMP,GTEMP CO2BD3 ! CO231,CO238,CDT31,CDT38 CO2BD2 ! C2D31,C2D38 CO2BD2 ! CO271,CO278,CDT71,CDT78 CO2BD4 ! C2D71,C2D78 CO2BD4 ! BETINW BDWIDE ! ! OUTPUTS: ! CO21,CO2NBL,CO2SP1,CO2SP2 TFCOM ! VAR1,VAR2,VAR3,VAR4,CNTVAL KDACOM ! QH2O,P,DELP,DELP2,T KDACOM ! ! CALLS: ! FST88 !----------------------------------------------------------------------- !----------------INPUT ARGUMENTS---------------------------------------- INTEGER, INTENT(IN), DIMENSION(:,:) :: KTOP,KBTM INTEGER, INTENT(IN), DIMENSION(:) :: NCLDS REAL, INTENT(IN), DIMENSION(:,:) :: EMCLD REAL, INTENT(IN), DIMENSION(:,:) :: PRES,TEMP,RH2O,QO3,CAMT REAL, INTENT(IN) :: RRVCO2 INTEGER, INTENT(IN),OPTIONAL, DIMENSION(:) :: LSFC REAL , INTENT(IN),OPTIONAL, DIMENSION(:) :: PSFC ! D I M E N S I O N ! & KTOP(IMAX,LP1),KBTM(IMAX,LP1),NCLDS(IMAX),EMCLD(IMAX,LP1), ! & TEMP(IMAX,LP1),PRES(IMAX,LP1),RH2O(IMAX,LMAX),QO3(IMAX,LMAX), ! & CAMT(IMAX,LP1),LSFC(IMAX),PSFC(IMAX) !----------------OUTPUT ARGUMENTS--------------------------------------- REAL, INTENT(OUT), DIMENSION(:,:) :: HEATRA REAL, INTENT(OUT), DIMENSION(:) :: GRNFLX,TOPFLX !----------------------------------------------------------------------- !----------------LOCAL ARRAY STORAGE------------------------------------ integer :: k, kp real :: diftt real, dimension(IMAX,LP1) :: & PRESS, CO2R1, DCO2D1, D2CD21, D2CD22, & CO2R2 , DCO2D2, TDAV, TSTDAV, VSUM3, TEXPSL, TLSQU real, dimension(IMAX,LMAX) :: & VV, CO2MR, CO2MD, CO2M2D, VSUM4 real, dimension(IMAX) :: & VSUM1, VSUM2, A1, A2 ! DIMENSION :: PRESS(IMAX,LP1) ! DIMENSION CO2R(IMAX,LP1,LP1),DIFT(IMAX,LP1,LP1) ! DIMENSION DIFT(IMAX,LP1,LP1) ! DIMENSION CO2R1 (IMAX,LP1), DCO2D1(IMAX,LP1) ! DIMENSION D2CD21(IMAX,LP1), D2CD22(IMAX,LP1) ! DIMENSION CO2R2 (IMAX,LP1), DCO2D2(IMAX,LP1) ! DIMENSION TDAV(IMAX,LP1),TSTDAV(IMAX,LP1), & ! VV(IMAX,LMAX),VSUM3(IMAX,LP1),VSUM1(IMAX),VSUM2(IMAX) ! DIMENSION CO2MR (IMAX,LMAX),CO2MD (IMAX,LMAX),CO2M2D(IMAX,LMAX) ! DIMENSION A1(IMAX),A2(IMAX) !----------------------------------------------------------------------- ! DIMENSION DIFTD(IMAX,LP1,LP1) ! DIMENSION DCO2DT(IMAX,LP1,LP1),D2CDT2(IMAX,LP1,LP1) ! DIMENSION TEXPSL(IMAX,LP1),TLSQU(IMAX,LP1) ! DIMENSION VSUM4(IMAX,LMAX) ! DIMENSION DIFT1D(IMAX,LP1M) !----------------------------------------------------------------------- ! EQUIVALENCE (DIFTD,CO2R) ! EQUIVALENCE (DCO2DT,D2CDT2,CO2R) ! EQUIVALENCE (VSUM3,TLSQU,TEXPSL) ! EQUIVALENCE (VV,VSUM4) ! EQUIVALENCE (DIFT1D,DIFT) !----------------------------------------------------------------------- !----------- check presence of optional arguments ---------------------- IF (( PRESENT(LSFC) .and. .not.PRESENT(PSFC)) .or. & (.not.PRESENT(LSFC) .and. PRESENT(PSFC))) THEN Call Error_Mesg ('LWRAD in LONGWAVE_MOD', & 'LSFC and PSFC must be used together.', FATAL) ENDIF !----------------------------------------------------------------------- !--- NOTE: convert input pressure (PRES) from MKS to CGS units DO K=1,LP1 DO I=1,IMAX PRESS(I,K)=10.*PRES(I,K) ENDDO ENDDO !****COMPUTE FLUX PRESSURES (P) AND DIFFERENCES (DELP2,DELP) DO K=2,LMAX DO I=1,IMAX P(I,K)=0.50*(PRESS(I,K-1)+PRESS(I,K)) ENDDO ENDDO IF (.not.PRESENT(LSFC)) THEN DO I=1,IMAX P(I,1)=0.0 P(I,LP1)=PRESS(I,LP1) ENDDO ELSE DO I=1,IMAX P(I,1)=0.0 P(I,LP1)=PRESS(I,LP1) !--- note: convert from MKS to CGS units P(I,LSFC(I)+1)=10.*PSFC(I) ENDDO ENDIF DO K=1,LMAX DO I=1,IMAX DELP2(I,K)=P(I,K+1)-P(I,K) ENDDO ENDDO DO K=1,LMAX DO I=1,IMAX DELP(I,K)=1.0/DELP2(I,K) ENDDO ENDDO !****COMPUTE FLUX LEVEL TEMPERATURES (T) AND CONTINUUM TEMPERATURE ! CORRECTIONS (TEXPSL) DO K=2,LMAX DO I=1,IMAX TTTT(I,K)=0.50*(TEMP(I,K-1)+TEMP(I,K)) ENDDO ENDDO IF (.not.PRESENT(LSFC)) THEN DO I=1,IMAX TTTT(I,1)=TEMP(I,1) TTTT(I,LP1)=TEMP(I,LP1) ENDDO ELSE DO I=1,IMAX TTTT(I,1)=TEMP(I,1) TTTT(I,LP1)=TEMP(I,LP1) TTTT(I,LSFC(I)+1)=TEMP(I,LP1) ENDDO ENDIF !****COMPUTE ARGUMENT FOR CONT.TEMP.COEFF. ! (THIS IS 1800.(1./TEMP-1./296.)) DO K=1,LP1 DO I=1,IMAX TEXPSL(I,K)=1800./TEMP(I,K)-6.081081081 ENDDO ENDDO !...THEN TAKE EXPONENTIAL DO K=1,LP1 DO I=1,IMAX TEXPSL(I,K)=EXP(TEXPSL(I,K)) ENDDO ENDDO !***COMPUTE WEIGHTED TEMPERATURE (TDAV) AND PRESSURE (TSTDAV) INTEGRALS !----------------------------------------------------------------------- !***COMPUTE OPTICAL PATHS FOR H2O AND O3, USING THE DIFFUSIVITY ! APPROXIMATION FOR THE ANGULAR INTEGRATION (1.66). OBTAIN THE ! UNWEIGHTED VALUES(VAR1,VAR3) AND THE WEIGHTED VALUES(VAR2,VAR4). ! THE QUANTITIES H3M4(.0003) AND H3M3(.003) APPEARING IN THE VAR2 AND ! VAR4 EXPRESSIONS ARE THE APPROXIMATE VOIGT CORRECTIONS FOR H2O AND ! O3,RESPECTIVELY. ! DO K=1,LMAX DO I=1,IMAX QH2O(I,K)=RH2O(I,K)*DIFFCTR ENDDO ENDDO !---VV IS THE LAYER-MEAN PRESSURE (IN ATM),WHICH IS NOT THE SAME AS ! THE LEVEL PRESSURE (PRESS) DO K=1,LMAX DO I=1,IMAX VV(I,K)=0.50*(P(I,K+1)+P(I,K))*P0INV ENDDO ENDDO DO K=1,LMAX DO I=1,IMAX VAR1(I,K)=DELP2(I,K)*QH2O(I,K)*GINV VAR3(I,K)=DELP2(I,K)*QO3(I,K)*DIFFCTR*GINV VAR2(I,K)=VAR1(I,K)*(VV(I,K)+0.0003) VAR4(I,K)=VAR3(I,K)*(VV(I,K)+0.003) ENDDO ENDDO ! COMPUTE OPTICAL PATH FOR THE H2O CONTINUUM, USING ROBERTS COEFFS. ! (BETINW),AND TEMP. CORRECTION (TEXPSL). THE DIFFUSIVITY FACTOR ! (WHICH CANCELS OUT IN THIS EXPRESSION) IS ASSUMED TO BE 1.66. THE ! USE OF THE DIFFUSIVITY FACTOR HAS BEEN SHOWN TO BE A SIGNIFICANT ! SOURCE OF ERROR IN THE CONTINUUM CALCS.,BUT THE TIME PENALTY OF ! AN ANGULAR INTEGRATION IS SEVERE. DO K=1,LMAX DO I=1,IMAX CNTVAL(I,K)=TEXPSL(I,K)*RH2O(I,K)*VAR2(I,K)*BETINW/ & (RH2O(I,K)+RATH2OMW) ENDDO ENDDO !***COMPUTE WEIGHTED TEMPERATURE (TDAV) AND PRESSURE (TSTDAV) INTEGRALS ! FOR USE IN OBTAINING TEMP. DIFFERENCE BET. SOUNDING AND STD. ! TEMP. SOUNDING (DIFT) DO I=1,IMAX TSTDAV(I,1)=0.0 TDAV(I,1)=0.0 ENDDO DO K=1,LP1 DO I=1,IMAX VSUM3(I,K)=TEMP(I,K)-STEMP(K) ENDDO ENDDO DO K=1,LMAX DO I=1,IMAX VSUM2(I)=GTEMP(K)*DELP2(I,K) VSUM1(I)=VSUM2(I)*VSUM3(I,K) TSTDAV(I,K+1)=TSTDAV(I,K)+VSUM2(I) TDAV(I,K+1)=TDAV(I,K)+VSUM1(I) ENDDO ENDDO !***COMPUTE DIFT !2 DO K=1,LP1 !2 DO KP=1,LP1 !2 DO I=1,IMAX ! DIFTD(I,KP,K)=TSTDAV(I,KP)-TSTDAV(I,K) !2 DIFT (I,KP,K)=TSTDAV(I,KP)-TSTDAV(I,K) !2 ENDDO !2 ENDDO ! ENDDO !2 DO K=1,LP1 !2 DO KP=1,LP1 !2 DO I=1,IMAX ! DIFTD(I,KP,K)=1.0/(DIFTD(I,KP,K)+1.E-20) !2 DIFT (I,KP,K)=1.0/(DIFT (I,KP,K)+1.E-20) !2 ENDDO !2 ENDDO !2 ENDDO !9 DO K=1,LP1 !9 DO KP=1,LP1 !9 DO I=1,IMAX ! DIFT(I,KP,K)=TDAV(I,KP)-TDAV(I,K) !d DIFT(I,KP,K)=(TDAV(I,KP)-TDAV(I,K))*DIFT(i,kp,k) !9 DIFT(I,KP,K)=(TDAV(I,KP)-TDAV(I,K))*(1.0/(tstdav(i,kp)- !9 & tstdav(i,k) + 1.E-20)) !9 ENDDO !9 ENDDO !9 ENDDO ! DO K=1,LP1 ! DO KP=1,LP1 ! DO I=1,IMAX ! DIFT(I,KP,K)=DIFT(I,KP,K)*DIFTD(I,KP,K) ! ENDDO ! ENDDO ! ENDDO DO K=1,LMAX DO I=1,IMAX VSUM4(I,K)=0.50*(VSUM3(I,K+1)+VSUM3(I,K)) ENDDO ENDDO !9 DO K=2,LP1 !9 DO I=1,IMAX !9 DIFT(I,K,K)=VSUM4(I,K-1) !9 ENDDO !9 ENDDO !****EVALUATE COEFFICIENTS FOR CO2 PRESSURE INTERPOLATION (A1,A2) ! ---- NOTE: PRESS converted from MKS to CGS DO I=1,IMAX A1(I)=(PRESS(I,LP1)-P0XZP8)/P0XZP2 A2(I)=(P0-PRESS(I,LP1))/P0XZP2 ENDDO !***PERFORM CO2 PRESSURE INTERPOLATION ON ALL INPUTTED TRANSMISSION ! FUNCTIONS AND TEMP. DERIVATIVES !---SUCCESSIVELY COMPUTING CO2R,DCO2DT AND D2CDT2 IS DONE TO SAVE ! STORAGE (AT A SLIGHT LOSS IN COMPUTATION TIME) DO K=1,LP1 DO KP=1,LP1 DO I=1,IMAX ! CO2R(I,KP,K)=A1(I)*CO251(KP,K)+A2(I)*CO258(KP,K) CO21(I,KP,K)=A1(I)*CO251(KP,K)+A2(I)*CO258(KP,K) ENDDO ENDDO ENDDO DO K=1,LP1 DO I=1,IMAX CO2R1(I,K)=A1(I)*CO231(K)+A2(I)*CO238(K) D2CD21(I,K)=A1(I)*C2D31(K)+A2(I)*C2D38(K) DCO2D1(I,K)=A1(I)*CDT31(K)+A2(I)*CDT38(K) CO2R2(I,K)=A1(I)*CO271(K)+A2(I)*CO278(K) D2CD22(I,K)=A1(I)*C2D71(K)+A2(I)*C2D78(K) DCO2D2(I,K)=A1(I)*CDT71(K)+A2(I)*CDT78(K) ENDDO ENDDO DO K=1,LP1 DO I=1,IMAX D2CD21(I,K)=1.E-3*D2CD21(I,K) DCO2D1(I,K)=1.E-2*DCO2D1(I,K) D2CD22(I,K)=1.E-3*D2CD22(I,K) DCO2D2(I,K)=1.E-2*DCO2D2(I,K) ENDDO ENDDO DO K=1,LMAX DO I=1,IMAX CO2MR(I,K)=A1(I)*CO2M51(K)+A2(I)*CO2M58(K) CO2MD(I,K)=A1(I)*CDTM51(K)+A2(I)*CDTM58(K) CO2M2D(I,K)=A1(I)*C2DM51(K)+A2(I)*C2DM58(K) ENDDO ENDDO DO K=1,LMAX DO I=1,IMAX CO2MD(I,K)=1.E-2*CO2MD(I,K) CO2M2D(I,K)=1.E-3*CO2M2D(I,K) ENDDO ENDDO !***COMPUTE CO2 TEMPERATURE INTERPOLATIONS FOR ALL BANDS,USING DIFT ! DO K=1,LP1 ! DO KP=1,LP1 ! DO I=1,IMAX ! CO21(I,KP,K)=CO2R(I,KP,K) ! ENDDO ! ENDDO ! ENDDO ! DO K=1,LP1 ! DO KP=1,LP1 ! DO I=1,IMAX ! DCO2DT(I,KP,K)=A1(I)*CDT51(KP,K)+A2(I)*CDT58(KP,K) ! ENDDO ! ENDDO ! ENDDO ! ! DO K=1,LP1 ! DO KP=1,LP1 ! DO I=1,IMAX ! DCO2DT(I,KP,K)=1.E-2*DCO2DT(I,KP,K) ! ENDDO ! ENDDO ! ENDDO DO K=1,LP1 DO KP=1,LP1 DO I=1,IMAX ! CO21(I,KP,K)=CO21(I,KP,K)+DIFT(I,KP,K)*DCO2DT(I,KP,K) !5 CO21(I,KP,K)=CO21(I,KP,K)+DIFT(I,KP,K)*(1.0E-2*(a1(I)* !5 & CDT51(KP,K)+A2(i)*CDT58(KP,K))) if ( (kp .ne. k) .or. k .eq. 1) then CO21(I,KP,K)=CO21(I,KP,K)+ & (Tdav(i,kp) - tdav(i,k))*(1.0/ (tstdav(i,kp) - & tstdav(i,k) + 1.E-20))* & (1.0E-2*(a1(I)* & CDT51(KP,K)+A2(i)*CDT58(KP,K))) else if ( (kp .eq. k) .and. (k .ge. 2) ) then CO21(I,KP,K)=CO21(I,KP,K)+ & vsum4(i,k-1)* & (1.0E-2*(a1(I)* & CDT51(KP,K)+A2(i)*CDT58(KP,K))) endif ENDDO ENDDO ENDDO ! DO K=1,LP1 ! DO KP=1,LP1 ! DO I=1,IMAX ! D2CDT2(I,KP,K)=A1(I)*C2D51(KP,K)+A2(I)*C2D58(KP,K) ! ENDDO ! ENDDO ! ENDDO ! DO K=1,LP1 ! DO KP=1,LP1 ! DO I=1,IMAX ! D2CDT2(I,KP,K)=1.E-3*D2CDT2(I,KP,K) ! ENDDO ! ENDDO ! ENDDO DO K=1,LP1 DO KP=1,LP1 DO I=1,IMAX ! CO21(I,KP,K)=CO21(I,KP,K)+0.50* ! DIFT(I,KP,K)* ! & DIFT(I,KP,K)*D2CDT2(I,KP,K) !7 CO21(I,KP,K)=CO21(I,KP,K)+0.50* !7 & DIFT(I,KP,K)* !7 & DIFT(I,KP,K)*(1.E-3*(A1(I)*C2D51(KP,K) + A2(I)* !7 & c2D58(KP,K))) if ( (kp .ne. k) .or. k .eq. 1) then DIFTT =(TDAV(I,KP)-TDAV(I,K))*(1.0/(tstdav(i,kp)- & tstdav(i,k) + 1.E-20)) ! CO21(I,KP,K)=CO21(I,KP,K)+0.50* CO21(I,KP,K)=CO21(I,KP,K)+0.50* & diftt*diftt* & !8 & DIFT(I,KP,K)* !8 & DIFT(I,KP,K)* ! & (0.5*((Tdav(i,kp) - tdav(i,k))*(1.0/ (tstdav(i,kp) - ! & tstdav(i,k) + 1.E-20)) * ! & (Tdav(i,kp) - tdav(i,k))* (1.0/ (tstdav(i,kp) - ! & tstdav(i,k) + 1.E-20)) ) )* ((1.E-3*(A1(I)*C2D51(KP,K) + A2(I)*c2D58(KP,K)))) ! & (1.E-3*(A1(I)*C2D51(KP,K) + A2(I)*c2D58(KP,K))) else if ( (kp .eq. k) .and. (k .ge. 2) ) then CO21(I,KP,K)=CO21(I,KP,K)+0.50* & vsum4(i,k-1)* & vsum4(i,k-1)* & (1.E-3*(A1(I)*C2D51(KP,K) + A2(I)*c2D58(KP,K))) endif ENDDO ENDDO ENDDO !***COMPUTE TRANSMISSION FCTNS USED IN SPA88 !---(IN THE 250 LOOP,DIFT REALLY SHOULD BE (I,1,K), BUT DIFT IS ! INVARIANT WITH RESPECT TO K,KP,AND SO (I,1,K)=(I,K,1)) DO K=1,LP1 DO I=1,IMAX !6 CO2SP1(I,K)=CO2R1(I,K)+DIFT(I,K,1)* !6 & (DCO2D1(I,K)+0.50*DIFT(I,K,1)*D2CD21(I,K)) !6 CO2SP2(I,K)=CO2R2(I,K)+DIFT(I,K,1)* !6 & (DCO2D2(I,K)+0.50*DIFT(I,K,1)*D2CD22(I,K)) CO2SP1(I,K)=CO2R1(I,K)+ & ((Tdav(i,k) - tdav(i,1))*(1.0/ (tstdav(i,k) - & tstdav(i,1) + 1.E-20))) * & (DCO2D1(I,K)+0.50* & ((Tdav(i,k ) - tdav(i,1))*(1.0/ (tstdav(i,k ) - & tstdav(i,1) + 1.E-20)))* & D2CD21(I,K)) CO2SP2(I,K)=CO2R2(I,K)+ & ((Tdav(i,k) - tdav(i,1)) *(1.0/ (tstdav(i,k ) - & tstdav(i,1) + 1.E-20)))* & (DCO2D2(I,K)+0.50* & ((Tdav(i,k) - tdav(i,1))*(1.0/ (tstdav(i,k ) - & tstdav(i,1) + 1.E-20)))* & D2CD22(I,K)) ENDDO ENDDO !--- WE AREN''T DOING NBL TFS ON THE 100 CM-1 BANDS . DO K=1,LMAX DO I=1,IMAX !7 CO2NBL(I,K)=CO2MR(I,K)+DIFT(I,K,K+1)*(CO2MD(I,K)+0.50* !7 & DIFT(I,K,K+1)*CO2M2D(I,K)) CO2NBL(I,K)=CO2MR(I,K)+ & ((Tdav(i,k ) - tdav(i,k+1))*(1.0/ (tstdav(i,k ) - & tstdav(i,k+1) + 1.E-20)))* & (CO2MD(I,K)+0.50* & ((Tdav(i,k ) - tdav(i,k+1))*(1.0/ (tstdav(i,k ) - & tstdav(i,k+1) + 1.E-20)))* & CO2M2D(I,K)) ENDDO ENDDO !***COMPUTE TEMP. COEFFICIENT BASED ON TTTT(K) (SEE REF.2) DO K=1,LP1 DO I=1,IMAX IF (TTTT(I,K) <= 250.) THEN TLSQU(I,K)=B0+(TTTT(I,K)-250.)* & (B1+(TTTT(I,K)-250.)* & (B2+B3*(TTTT(I,K)-250.))) ELSE TLSQU(I,K)=B0 ENDIF ENDDO ENDDO !***APPLY TO ALL CO2 TFS DO K=1,LP1 DO KP=1,LP1 DO I=1,IMAX CO21(I,KP,K)=CO21(I,KP,K)*(1.0-TLSQU(I,KP))+TLSQU(I,KP) ENDDO ENDDO ENDDO DO K=1,LP1 DO I=1,IMAX CO2SP1(I,K)=CO2SP1(I,K)*(1.0-TLSQU(I,1))+TLSQU(I,1) CO2SP2(I,K)=CO2SP2(I,K)*(1.0-TLSQU(I,1))+TLSQU(I,1) ENDDO ENDDO DO K=1,LMAX DO I=1,IMAX CO2NBL(I,K)=CO2NBL(I,K)*(1.0-TLSQU(I,K))+TLSQU(I,K) ENDDO ENDDO !----------------------------------------------------------------------- IF (Present(LSFC)) THEN CALL FST88 (KTOP,KBTM,NCLDS,EMCLD,PRESS,TEMP,CAMT,RRVCO2, & HEATRA,GRNFLX,TOPFLX, LSFC) ELSE CALL FST88 (KTOP,KBTM,NCLDS,EMCLD,PRESS,TEMP,CAMT,RRVCO2, & HEATRA,GRNFLX,TOPFLX) ENDIF !----------------------------------------------------------------------- END SUBROUTINE LWRAD2 !####################################################################### !####################################################################### SUBROUTINE FST88 (KTOP,KBTM,NCLDS,EMCLD,PRESS,TEMP,CAMT,RRVCO2, & HEATRA,GRNFLX,TOPFLX, LSFC) !----------------------------------------------------------------------- ! INPUT PARAMETERS ! ---------------- ! ! KTOP = INDEX OF (DATA LEVEL) PRESSURE OF CLOUD TOP,USED ! IN THE LONGWAVE PROGRAM ! KBTM = INDEX OF (DATA LEVEL) PRESSURE OF CLOUD BOTTOM, ! USED IN THE LONGWAVE PROGRAM ! NCLDS = NO. CLOUDS AT EACH GRID PT. ! EMCLD = CLOUD EMISSIVITY. SET TO ONE BY DEFAULT, BUT MAY ! BE MODIFIED FOR USE IN LONGWAVE PROGRAM. ! PRESS = PRESSURE (CGS UNITS) AT DATA LEVELS OF MODEL ! TEMP = TEMPERATURE (K) AT DATA LEVELS OF MODEL ! CAMT = CLOUD AMOUNTS OF CLOUDS (THEIR LOCATIONS ARE ! SPECIFIED IN THE KTOP/KBTM INDICES) ! RRVCO2 = THE VOLUME MIXING RATIO OF CO2 (SCALAR) ! (used in NLTE, passed by argument list) ! ! LSFC = Vertical index of the lowest model level, ! dimensioned by IMAX. (For eta coordinate) ! !----------------------------------------------------------------------- ! OUTPUT PARAMETERS ! ----------------- ! ! HEATRA = HEATING RATE AT DATA LEVELS (K/DAY) ! GRNFLX = NET LONGWAVE FLUX AT THE GROUND (CGS UNITS) ! TOPFLX = NET LONGWAVE FLUX AT THE TOP (CGS UNITS) ! ! !----------------------------------------------------------------------- ! ***************************************************************** ! SUBROUTINE FST88 IS THE MAIN COMPUTATION MODULE OF THE ! LONG-WAVE RADIATION CODE. IN IT ALL "EMISSIVITY" CALCULATIONS, ! INCLUDING CALLS TO TABLE LOOKUP SUBROUTINES. ALSO,AFTER CALLING ! SUBROUTINE "SPA88", FINAL COMBINED HEATING RATES AND GROUND ! FLUX ARE OBTAINED. ! ***************************************************************** !----------------------------------------------------------------------- ! OTHER INPUTS: ! ------------- ! BETINW,BETAWD,AB15WD BDWIDE ! BETAD,BO3RND,AO3RND BANDTA ! QH2O,P,DELP2,DELP,T,VAR1,VAR2, KDACOM ! VAR3,VAR4,CNTVAL KDACOM ! BCLDS RDBITS ! INDX2,KMAXV,SOURCE,DSRCE TABCOM ! SKC1R,SKC3R,KMAXVM,NREP1,NREP2 TABCOM ! CO2NBL,CO2SP,CO21 TFCOM ! ! OUTPUTS: ! FLX1E1 RDFLUX ! ! CALLED BY : LWR88 ! ! CALLS : CLO88,E1E288,E3V88,SPA88,NLTE ! ! ! PASSED VARIABLES: ! ! IN E3V88: ! EMD = E3 FUNCTION FOR H2O LINES (0-560,1200-2200 CM-1) ! COMPUTED IN E3V88 ! TPL = TEMPERATURE INPUT FOR E3 CALCULATION IN E3V88 ! EMPL = H2O AMOUNT,INPUT FOR E3 CALCULATION IN E3V88 ! ! IN E1E288: ! E1CTS1 = E1 FUNCTION FOR THE (I+1)TH LEVEL USING THE ! TEMPERATURE OF THE ITH DATA LEVEL,COMPUTED OVER ! THE FREQUENCY RANGE 0-560,1200-2200 CM-1. (E1CTS1- ! E1CTW1) IS USED IN OBTAINING THE FLUX AT THE TOP ! IN THE 0-160,1200-2200 CM-1 RANGE (FLX1E1). ! E1CTS2 = E1 FUNCTION FOR THE ITH LEVEL, USING THE TEMP. OF ! THE ITH DATA LEVEL,COMPUTED OVER THE FREQUENCY RANGE ! 0-560,1200-2200 CM-1. (E1CTS2-E1CTW2) IS ALSO USED ! IN OBTAINING THE FLUX AT THE TOP IN THE 0-160,. ! 1200-2200 CM-1 RANGE. ! E1FLX = E1 FCTN. FOR THE ITH LEVEL,USING THE TEMPERATURE AT ! THE TOP OF THE ATMOSPHERE. COMPUTED OVER THE FREQ. ! RANGE 0-560,1200-2200 CM-1. USED FOR Q(APPROX) TERM. ! (IN COMMON BLOCK TFCOM) ! E1CTW1 = LIKE E1CTS1,BUT COMPUTED OVER THE 160-560 CM-1 RANGE ! AND USED FOR Q(APPROX,CTS) CALCULATION ! E1CTW2 = LIKE E1CTS2,BUT COMPUTED OVER THE 160-560 CM-1 RANGE ! AND USED FOR Q(APPROX,CTS) CALCULATION ! FXO = TEMPERATURE INDEX USED FOR E1 FUNCTION AND ALSO ! USED FOR SOURCE FUNCTION CALC. IN FST88. ! DT = TEMP. DIFF.BETWEEN MODEL TEMPS. AND TEMPS. AT ! TABULAR VALUES OF E1 AND SOURCE FCTNS. USED IN ! FST88 AND IN E1 FUNCTION CALC. !----------------------------------------------------------------------- !#include "hcon.h" !----------------------------------------------------------------------- !#include "bandta.h" !----------------------------------------------------------------------- !#include "bdwide.h" !----------------------------------------------------------------------- !#include "bdcomb.h" !----------------------------------------------------------------------- !----------------INPUT ARGUMENTS---------------------------------------- INTEGER, INTENT(IN), DIMENSION(:,:) :: KTOP,KBTM INTEGER, INTENT(IN), DIMENSION(:) :: NCLDS REAL, INTENT(IN), DIMENSION(:,:) :: EMCLD,PRESS,TEMP,CAMT REAL, INTENT(IN) :: RRVCO2 INTEGER, INTENT(IN), DIMENSION(:), OPTIONAL :: LSFC ! D I M E N S I O N ! & KTOP(IMAX,LP1),KBTM(IMAX,LP1),NCLDS(IMAX),EMCLD(IMAX,LP1), ! & TEMP(IMAX,LP1),PRESS(IMAX,LP1),CAMT(IMAX,LP1),LSFC(IMAX) !----------------OUTPUT ARGUMENTS--------------------------------------- REAL, INTENT(OUT), DIMENSION(:,:) :: HEATRA REAL, INTENT(OUT), DIMENSION(:) :: GRNFLX,TOPFLX !----------------OUTPUT ARGUMENTS--------------------------------------- ! DIMENSION HEATRA(IMAX,LMAX),GRNFLX(IMAX),TOPFLX(IMAX) !----------------------------------------------------------------------- !----------------LOCAL ARRAY STORAGE------------------------------------ integer :: i, k, n, klen, kmaxs, kk, kp, icnt, nc real :: hm5 integer :: IXO(IMAX,LP1),ITOP(IMAX),IBOT(IMAX) real :: VTMP1(IMAX,LP1M),VTMP2(IMAX,LP1V),VTMP3(IMAX,LP1), & C(IMAX,LLP1),ALP(IMAX,LLP1),DSORC(IMAX,LP1), & TOTPHI(IMAX,LP1),TOTO3(IMAX,LP1),TPHIO3(IMAX,LP1), & RLOG(IMAX,LMAX),DELPTC(IMAX),PTOP(IMAX),PBOT(IMAX), & FTOP(IMAX),FBOT(IMAX) !----------------------------------------------------------------------- ! ------DIMENSION OF VARIABLES EQUIVALENCED TO THOSE ABOVE------ real :: TVAL(IMAX,LP1),VSUM1(IMAX,LP1),HEATEM(IMAX,LP1) real :: EMXX(IMAX,LMAX) real :: CSUB(IMAX,LLP1) real :: FLX(IMAX,LP1) real :: OSS(IMAX,LP1),CSS(IMAX,LP1),SS2(IMAX,LP1),TC(IMAX,LP1), & DTC(IMAX,LP1) real :: ALPSQ1(IMAX,LP1),ALPSQ2(IMAX,LP1) real :: DELPR1(IMAX,LP1),DELPR2(IMAX,LP1) real :: FLXNET(IMAX,LP1),FLXTHK(IMAX,LP1) real :: Z1(IMAX,LP1),CEVAL(IMAX,LP1) real :: TOTEVV(IMAX,LP1) real :: AVMO3(IMAX,LP1M),AVPHO3(IMAX,LP1V),FAC1(IMAX,LP1V) real :: FAC2(IMAX,LP1M),AVVO2(IMAX,LP1M),AVEPHJ(IMAX,LP1V) !1 DIMENSION OVER(IMAX,LP1,LP1) ! DIMENSION OVER1D(IMAX,LP1M) !1 DIMENSION EMISST(IMAX,LP1,LP1) !---DIMENSION OF VARIABLES EQUIVALENCED TO THOSE IN OTHER COMMON BLKS-- ! DIMENSION TO31D(IMAX,LP1M),EMI21D(IMAX,LP1M) ! DIMENSION CO21D(IMAX,LP1M),EMIS1D(IMAX,LP1M) ! DIMENSION AVEP1D(IMAX,LP1M),AVEP1(IMAX*LP1M) !---DIMENSION OF VARIABLES PASSED TO OTHER SUBROUTINES--- real :: E1CTS1(IMAX,LP1),E1CTS2(IMAX,LMAX) real :: E1CTW1(IMAX,LP1),E1CTW2(IMAX,LMAX) real :: FXO(IMAX,LP1),DT(IMAX,LP1) real :: EMD(IMAX,LLP1),TPL(IMAX,LLP1),EMPL(IMAX,LLP1) !---EMX1 IS A LOCAL VARIABLE USED AS INPUT AND OUTPUT TO E1E288-- real :: EMX1(IMAX) ! EQUIVALENCE (VTMP3,TVAL,VSUM1,EMXX,TC,HEATEM) ! EQUIVALENCE (VTMP2,VTMP21,AVPHO3,FAC1,AVEPHJ) ! EQUIVALENCE (VTMP1,AVMO3,FAC2,AVVO2,OVER,EMISST) ! EQUIVALENCE (DSORC,TOTEVV,DELPR1,OSS,SUM,FLXNET) ! EQUIVALENCE (TOTPHI,DELPR2,CSS,FLX,Z1) ! EQUIVALENCE (TOTO3,ALPSQ1,DTC,FLXTHK) ! EQUIVALENCE (TPHIO3,ALPSQ2,SS2,CEVAL) ! EQUIVALENCE (AVEPHI,AVEP1D,AVEP1),(EMI21D,EMISS2) ! EQUIVALENCE (OVER1D,OVER),(TO31D,TO3),(EMIS1D,EMISS) ! EQUIVALENCE (CO21D,CO21) ! EQUIVALENCE (ALP,CSUB) !----------------------------------------------------------------------- ! FIRST SECTION IS TABLE LOOKUP FOR SOURCE FUNCTION AND ! DERIVATIVE (B AND DB/DT).ALSO,THE NLTE CO2 SOURCE FUNCTION ! IS OBTAINED DO K=1,LP1 DO I=1,IMAX FXO(I,K)=AINT((TEMP(I,K)-100.)*0.10+1.0) ENDDO ENDDO DO K=1,LP1 DO I=1,IMAX TVAL(I,K)=90.+10.*FXO(I,K) ENDDO ENDDO DO K=1,LP1 DO I=1,IMAX DT(I,K)=TEMP(I,K)-TVAL(I,K) ENDDO ENDDO DO K=1,LP1 DO I=1,IMAX IXO(I,K)=FXO(I,K) SS1(I,K)=0.0 ENDDO ENDDO !----------------------------------------------------------------------- DO 110 N=1,NBLY DO K=1,LP1 DO I=1,IMAX SORC(I,K,N)=SOURCE(IXO(I,K),N) DSORC(I,K)=DSRCE(IXO(I,K),N) ENDDO ENDDO DO K=1,LP1 DO I=1,IMAX SORC(I,K,N)=SORC(I,K,N)+DT(I,K)*DSORC(I,K) ENDDO ENDDO IF (N == 11 .or. N == 12 .or. N == 14) THEN DO K=1,LP1 DO I=1,IMAX SS1(I,K)=SS1(I,K)+SORC(I,K,N) ENDDO ENDDO ENDIF 110 CONTINUE !----------------------------------------------------------------------- DO K=1,LP1 DO I=1,IMAX CSOUR1(I,K)=SORC(I,K,9) ENDDO ENDDO DO K=1,LP1 DO I=1,IMAX CSOUR2(I,K)=SORC(I,K,10) ENDDO ENDDO DO K=1,LP1 DO I=1,IMAX OSOUR(I,K)=SORC(I,K,13) ENDDO ENDDO !----------------------------------------------------------------------- ! ! THE FOLLOWING SUBROUTINE OBTAINS NLTE SOURCE FUNCTION FOR CO2 ! ! CALL NLTE (PRESS,RRVCO2) ! DO 141 K=1,INLTE ! DO 141 I=1,IMAX ! CSOUR1(I,K)=SORC(I,K,9) ! CSOUR2(I,K)=SORC(I,K,10) !141 CONTINUE ! !----------------------------------------------------------------------- DO K=1,LP1 DO I=1,IMAX CSOUR(I,K)=CSOUR1(I,K)+CSOUR2(I,K) ENDDO ENDDO !----------------------------------------------------------------------- ! SECOND SECTION PRODUCES 4 MATRICES FOR LATER USE IN PROGRAM: ! 1) AVEPHI(I,J) IS THE SCALED WATER MASS FOR USE IN TABLES (THIS ! IS U(P,P'') IN REF.(4); ! 2) OVER (I,J) IS THE WATER TRANSMISSION FUNCTION (USING ! "EMISSIVITY" APPROXIMATION) IN THE 560-800 CM-1 BAND; ! 3) TO3(I,J) IS THE EXACT OZONE TRANSMISSION FUNCTION (USING ! PARAMETERS FOR THE 990-1070 CM-1 BAND FROM THE 1982 AFGL CATALOG) ! 4)EMISS2(I,J) IS THE EMISSIVITY H20 TRANSMISSION FUNCTION DUE ! TO THE 10 UM CONTINUUM,TREATED AS ONE BAND. DO I=1,IMAX TOTPHI(I,1)=0.0 TOTO3(I,1)=0.0 TPHIO3(I,1)=0.0 TOTVO2(I,1)=0.0 ENDDO DO K=2,LP1 DO I=1,IMAX TOTPHI(I,K)=TOTPHI(I,K-1)+VAR2(I,K-1) ENDDO ENDDO DO K=2,LP1 DO I=1,IMAX TOTO3(I,K)=TOTO3(I,K-1)+VAR3(I,K-1) ENDDO ENDDO DO K=2,LP1 DO I=1,IMAX TPHIO3(I,K)=TPHIO3(I,K-1)+VAR4(I,K-1) ENDDO ENDDO DO K=2,LP1 DO I=1,IMAX TOTVO2(I,K)=TOTVO2(I,K-1)+CNTVAL(I,K-1) ENDDO ENDDO ! THE CALCULATIONAL PROCEDURE USED HERE IS: ! 1) FROM THE 1-D VECTORS (TOTPHI,TOTO3,ETC) CONSTRUCT A 1-D ! ARRAYS INCLUDING ALL INFORMATION FOR THE UPPER (I>J) TRIANGLE ! OF THE SYMMETRIC MATRICES(LP1,LP1); ! 2)PERFORM COMPUTATIONS ON THESE ARRAYS TO OBTAIN TO3,OVER, ! AVEPHJ: UPPER TRIANGLE OF THE SYMMETRIC MATRICES IN 1-D FORM ! 3) LOAD NUMBERS INTO THE UPPER TRIANGLE OF THE 2-D ! MATRICES TO3,OVER,AVEPHI,USING THE CDC Q8MERG FUNCTION ! 4) FILL UP THE LOWER TRIANGLE,USING THE CDC "SCATTER" FUNCTION ! THE DIAGRAM BELOW ILLUSTRATES THE RELATIONSHIP BETWEEN THE 1-D ! ARRAY AND THE 2-D SYMMETRIC MATRIX FOR A 4X4 MATRIX. ! ! I ! 1 2 3 4 ! -------------------------- ! 1 1 2 3 THE NOS. ARE THE ! J 2 4 5 POSITIONS IN THE ! 3 6 1-D ARRAY ! 4 ! !***IN ORDER TO USE A MINIMUM OF STORAGE SPACE, THE CALCULATIONS HAVE ! BEEN ORDERED AS FOLLOWS: ! ! --------------------------------------- ! STAGE 1...COMPUTE O3 TRANSMISSION FCTNS ! --------------------------------------- DO K=1,LMAX KLEN=LP1-K KMAXS=KMAXV(K) DO KK=1,KLEN DO I=1,IMAX AVMO3(I,KMAXS+KK-1)=TOTO3(I,KK+K)-TOTO3(I,K) AVPHO3(I,KMAXS+KK-1)=TPHIO3(I,KK+K)-TPHIO3(I,K) ENDDO ENDDO ENDDO !---SOME OF THE CALCS. BELOW ARE "STRIP-MINED" SO THAT VECTOR CALCS. ! ARE RESTRICTED TO A LENGTH OF 65535. THIS IS A CYBER 205 LIMIT. ! (THE COMPILER HANDLES STRIP-MINING OF FORTRAN 77 CODE). DO KK=1,KMAXVM DO I=1,IMAX FAC1(I,KK)=BO3RND(2)*AVPHO3(I,KK)/AVMO3(I,KK) ENDDO ENDDO DO KK=1,KMAXVM DO I=1,IMAX FAC2(I,KK)=4.0*AO3RND(2)*AVMO3(I,KK) ENDDO ENDDO DO KK=1,KMAXVM DO I=1,IMAX VTMP1(I,KK)=1.0+FAC2(I,KK)/FAC1(I,KK) ENDDO ENDDO DO KK=1,KMAXVM DO I=1,IMAX VTMP1(I,KK)=SQRT(VTMP1(I,KK)) ENDDO ENDDO DO KK=1,KMAXVM DO I=1,IMAX VTMP2(I,KK)=0.50*(FAC1(I,KK)*(VTMP1(I,KK)-1.0)) ENDDO ENDDO DO K=1,LMAX KLEN=LP1-K KMAXS=KMAXV(K) DO KK=1,KLEN DO I=1,IMAX AVVO2(I,KMAXS+KK-1)=TOTVO2(I,KK+K)-TOTVO2(I,K) ENDDO ENDDO ENDDO DO K=1,LMAX DO I=1,IMAX TO3SPC(I,K)=VTMP2(I,K) ENDDO ENDDO DO KK=1,KMAXVM DO I=1,IMAX VTMP2(I,KK)=VTMP2(I,KK)+SKO3R*AVVO2(I,KK) ENDDO ENDDO DO KK=1,KMAXVM DO I=1,IMAX VTMP2(I,KK)=-1.0*VTMP2(I,KK) ENDDO ENDDO DO KK=1,KMAXVM DO I=1,IMAX VTMP2(I,KK)=EXP(VTMP2(I,KK)) ENDDO ENDDO DO K=1,LMAX KMAXS=KMAXV(K) DO KP=K+1,LP1 DO I=1,IMAX TO3(I,KP,K)=VTMP2(I,KMAXS-(K+1)+KP) ENDDO ENDDO ENDDO DO K=1,LP1 DO I=1,IMAX TO3(I,K,K)=1.0 ENDDO ENDDO DO K=1,KMAXVM DO I=1,IMAX TO3(I,INDX1(K),INDX2(K))=VTMP2(I,K) ENDDO ENDDO ! ! ------------------------- ! STAGE 2... COMPUTE AVEPHI ! ------------------------- ! DO 401 K=1,LMAX KLEN=LP1-K KMAXS=KMAXV(K) DO KK=1,KLEN DO I=1,IMAX AVEPHJ(I,KMAXS+KK-1)=TOTPHI(I,KK+K)-TOTPHI(I,K) ENDDO ENDDO DO KP=K+1,LP1 DO I=1,IMAX AVEPHI(I,KP,K)=AVEPHJ(I,KMAXS-(K+1)+KP) ENDDO ENDDO 401 CONTINUE DO K=1,LP1 DO I=1,IMAX AVEPHI(I,K,K)=1.01E-16 ENDDO ENDDO DO K=1,KMAXVM DO I=1,IMAX AVEPHI(I,INDX1(K),INDX2(K))=AVEPHJ(I,K) ENDDO ENDDO ! ---------------------- ! STAGE 3...COMPUTE OVER ! ---------------------- DO KK=1,KMAXVM DO I=1,IMAX VTMP2(I,KK)=AB15WD*AVEPHJ(I,KK) ENDDO ENDDO DO KK=1,KMAXVM DO I=1,IMAX VTMP2(I,KK)=SQRT(VTMP2(I,KK)) ENDDO ENDDO DO KK=1,KMAXVM DO I=1,IMAX VTMP2(I,KK)=VTMP2(I,KK)+SKC1R*AVVO2(I,KK) ENDDO ENDDO DO KK=1,KMAXVM DO I=1,IMAX VTMP2(I,KK)=-1.0*VTMP2(I,KK) ENDDO ENDDO DO KK=1,KMAXVM DO I=1,IMAX VTMP2(I,KK)=EXP(VTMP2(I,KK)) ENDDO ENDDO ! DO K=1,LMAX ! KMAXS=KMAXV(K) ! DO KP=K+1,LP1 ! DO I=1,IMAX ! OVER(I,KP,K)=VTMP2(I,KMAXS-(K+1)+KP) ! ENDDO ! ENDDO ! ENDDO ! DO K=1,LP1 ! DO I=1,IMAX ! OVER(I,K,K)=1.0 ! ENDDO ! ENDDO ! DO K=1,KMAXVM ! DO I=1,IMAX ! OVER(I,INDX1(K),INDX2(K))=VTMP2(I,K) ! ENDDO ! ENDDO ! ------------------------ ! STAGE 4...COMPUTE EMISS2 ! ------------------------ DO K=1,LP1 DO I=1,IMAX VTMP3(I,K)=TOTVO2(I,K) ENDDO ENDDO DO K=1,LP1 DO I=1,IMAX VTMP3(I,K)=-1.0*VTMP3(I,K) ENDDO ENDDO DO K=1,LP1 DO I=1,IMAX VTMP3(I,K)=EXP(VTMP3(I,K)) ENDDO ENDDO ! COMPUTE 1-BAND CONTINUUM TRANSMISSIVITIES (EMISS2) DO K=1,LP1 DO I=1,IMAX TOTEVV(I,K)=1.0/VTMP3(I,K) ENDDO ENDDO DO K=1,LP1 DO KP=1,K DO I=1,IMAX EMISS2(I,KP,K)=VTMP3(I,K)*TOTEVV(I,KP) ENDDO ENDDO ENDDO DO K=1,LMAX DO KP=K+1,LP1 DO I=1,IMAX EMISS2(I,KP,K)=VTMP3(I,KP)*TOTEVV(I,K) ENDDO ENDDO ENDDO !----------------------------------------------------------------------- ! THE THIRD SECTION CALCULATES BOUNDARY LAYER AND NEARBY LAYER ! CORRECTIONS TO THE TRANSMISSION FUNCTIONS OBTAINED ABOVE. METHODS ! ARE GIVEN IN REF. (4). !----------------------------------------------------------------------- ! THE FOLLOWING RATIOS ARE USED IN VARIOUS NBL CALCULATIONS: DO K=2,LMAX DO I=1,IMAX DELPR1(I,K)=DELP(I,K)*(PRESS(I,K)-P(I,K)) ENDDO ENDDO DO K=2,LP1 DO I=1,IMAX DELPR2(I,K)=DELP(I,K-1)*(P(I,K)-PRESS(I,K-1)) ENDDO ENDDO ! THE FIRST COMPUTATION IS FOR THE 15 UM BAND,WITH THE H2O ! TRANS. FCTN (OVER) COMBINED WITH THE CO2 FCTN (CO21). ! ! COMBINE CO21,OVER INTO CO21; BEFORE MAKING NBL CORRECTIONS, ! LOAD (1,K) VALUES FOR USE IN EXACT CTS CALCULATIONS IN SPA88. DO K=1,LMAX KMAXS=KMAXV(K) DO KP=K+1,LP1 DO I=1,IMAX CO21(I,KP,K)=CO21(I,KP,K)*VTMP2(I,KMAXS-(K+1)+KP) ENDDO ENDDO ENDDO ! DO K=1,LP1 ! DO I=1,IMAX ! OVER(I,K,K)=1.0 ! ENDDO ! ENDDO DO K=1,KMAXVM DO I=1,IMAX ! CO21(I,KP,K)=CO21(I,KP,K)* CO21(I,indx1(k), indx2(k))=CO21(I,indx1(k), indx2(k))* & VTMP2(I,K) ENDDO ENDDO ! DO K=1,LP1 ! DO KP=1,LP1 ! DO I=1,IMAX ! CO21(I,KP,K)=CO21(I,KP,K)*OVER(I,KP,K) ! ENDDO ! ENDDO ! ENDDO DO K=1,LP1 DO I=1,IMAX CO2SP(I,K)=CO21(I,1,K) ENDDO ENDDO ! PERFORM NBL COMPUTATIONS FOR THE CO2 15 UM BAND DO K=2,LMAX DO I=1,IMAX ALPSQ1(I,K)=SQRT(DELPR1(I,K)) ENDDO ENDDO DO K=2,LP1 DO I=1,IMAX ALPSQ2(I,K)=SQRT(DELPR2(I,K)) ENDDO ENDDO ! DO K=1,LMAX ! KMAXS=KMAXV(K) ! DO KP=K+1,LP1 ! DO I=1,IMAX ! OVER(I,KP,K)=VTMP2(I,KMAXS-(K+1)+KP) ! ENDDO ! ENDDO ! ENDDO ! DO K=1,LP1 ! DO I=1,IMAX ! OVER(I,K,K)=1.0 ! ENDDO ! ENDDO ! DO K=1,KMAXVM ! DO I=1,IMAX ! OVER(I,INDX1(K),INDX2(K))=VTMP2(I,K) ! ENDDO ! ENDDO ! DO K=1,LMAX ! DO I=1,IMAX ! RLOG(I,K)=OVER(I,K,K+1) ! ENDDO ! ENDDO do k=1,LMAX do i=1,IMAX do kk=1,KMAXVM if (indx1(kk) .eq. k .and. indx2(kk) .eq. (k+1) ) then rlog(i,k) = vtmp2(i,kk) cycle endif end do end do end do DO K=1,LMAX DO I=1,IMAX RLOG(I,K)=RLOG(I,K)*CO2NBL(I,K) ENDDO ENDDO DO K=1,LMAX DO I=1,IMAX RLOG(I,K)=LOG(RLOG(I,K)) ENDDO ENDDO DO K=2,LMAX DO I=1,IMAX ALP(I,LM1+K)=-ALPSQ1(I,K)*RLOG(I,K) ENDDO ENDDO DO K=2,LP1 DO I=1,IMAX ALP(I,K-1)=-ALPSQ2(I,K)*RLOG(I,K-1) ENDDO ENDDO DO I=1,IMAX ALP(I,LL)=-RLOG(I,LMAX) ALP(I,LLP1)=-RLOG(I,LMAX)* & SQRT(DELP(I,LMAX)*(P(I,LP1)-PRESS(I,LM1))) ENDDO !***THE STATEMENT FUNCTION SF IN PREV. VERSIONS IS NOW EXPLICITLY ! EVALUATED. DO K=1,LLP1 DO I=1,IMAX C(I,K)=ALP(I,K)*(-0.66667+ALP(I,K)*(0.25+ALP(I,K)*(-0.066667))) ENDDO ENDDO DO I=1,IMAX CO21(I,LP1,LP1)=1.0+C(I,LMAX) CO21(I,LP1,LMAX)=1.0+(DELP2(I,LMAX)*C(I,LL)-(PRESS(I,LMAX)- & P(I,LMAX))*C(I,LLM1))/(P(I,LP1)-PRESS(I,LMAX)) CO21(I,LMAX,LP1)=1.0+((P(I,LP1)-PRESS(I,LM1))*C(I,LLP1)- & (P(I,LP1)-PRESS(I,LMAX))*C(I,LMAX))/ & (PRESS(I,LMAX)-PRESS(I,LM1)) ENDDO !*** THE K INDICES IN THE FOLLOWING LOOP RUN FROM 21 TO 40 IN THE ! L40 SKYHI CODE VERSION DO K=2,LMAX DO I=1,IMAX CEVAL(I,K-1)=1.0+0.50*(C(I,LM1+K)+C(I,K-1)) ENDDO ENDDO DO K=2,LMAX DO I=1,IMAX CO21(I,K,K)=CEVAL(I,K-1) ENDDO ENDDO ! COMPUTE NEARBY-LAYER TRANSMISSIVITIES FOR THE O3 BAND AND FOR THE ! ONE-BAND CONTINUUM BAND (EMISS2 AND TO3). THE SF2 FUNCTION IS ! USED. THE METHOD IS THE SAME AS DESCRIBED FOR CO2 IN REF (4). DO K=2,LMAX DO I=1,IMAX CSUB(I,K)=CNTVAL(I,K)*DELPR1(I,K) ENDDO ENDDO DO K=2,LMAX DO I=1,IMAX CSUB(I,K+LM1)=CNTVAL(I,K-1)*DELPR2(I,K) ENDDO ENDDO !---THE SF2 FUNCTION IN PREV. VERSIONS IS NOW EXPLICITLY EVALUATED HM5=-0.50 DO K=2,LLM1 DO I=1,IMAX ! C(I,K)=CSUB(I,K)*(HM5+CSUB(I,K)* ! & (0.166666-CSUB(I,K)*0.0416666)) C(I,K)= 0.166666-CSUB(I,K)*0.0416666 C(I,K)= CSUB(I,K)*(CSUB(I,K)*C(I,K)-0.50) ENDDO ENDDO DO I=1,IMAX EMISS2(I,LP1,LP1)=1.0+C(I,LLM1) ENDDO DO K=2,LMAX DO I=1,IMAX CEVAL(I,K-1)=1.0+0.50*(C(I,K)+C(I,LM1+K)) ENDDO ENDDO DO K=2,LMAX DO I=1,IMAX EMISS2(I,K,K)=CEVAL(I,K-1) ENDDO ENDDO DO K=2,LLM1 DO I=1,IMAX CSUB(I,K)=CSUB(I,K)*SKO3R ENDDO ENDDO ! !---THE SF2 FUNCTION IN PREV. VERSIONS IS NOW EXPLICITLY EVALUATED DO K=2,LLM1 DO I=1,IMAX C(I,K)=CSUB(I,K)*(-0.50+CSUB(I,K)* & (0.166666-CSUB(I,K)*0.0416666)) ENDDO ENDDO DO I=1,IMAX TO3(I,LP1,LP1)=1.0+C(I,LLM1) ENDDO DO K=2,LMAX DO I=1,IMAX CEVAL(I,K-1)=1.0+0.50*(C(I,K)+C(I,LM1+K)) ENDDO ENDDO DO K=2,LMAX DO I=1,IMAX TO3(I,K,K)=CEVAL(I,K-1) ENDDO ENDDO !----------------------------------------------------------------------- ! FOURTH SECTION OBTAINS WATER TRANSMISSION FUNCTIONS ! USED IN Q(APPROX) CALCULATIONS AND ALSO MAKES NBL CORRECTIONS: ! 1) EMISS (I,J) IS THE TRANSMISSION FUNCTION MATRIX OBTAINED ! BY CALLING SUBROUTINE E1E288; ! 2) "NEARBY LAYER" CORRECTIONS (EMISS(I,I)) ARE OBTAINED ! USING SUBROUTINE E3V88; ! 3) SPECIAL VALUES AT THE SURFACE (EMISS(LMAX,LP1),EMISS(LP1,L), ! EMISS(LP1,LP1)) ARE CALCULATED. ! ! EMXX,AV1,AND EMPL ARE COMPUTED BEFORE AVEPHI IS MODIFIED ! ! OBTAIN ARGUMENTS FOR E1E288 AND E3V88: !----------------------------------------------------------------------- ! DO I=1,IMAX EMX1(I)=QH2O(I,LMAX)*PRESS(I,LMAX)*(PRESS(I,LMAX)-P(I,LMAX))* & GP0INV ENDDO DO K=1,LM1 DO I=1,IMAX EMXX(I,K)=AVEPHI(I,K,LMAX)+EMX1(I) ENDDO ENDDO DO I=1,IMAX EMPL(I,1)=AVEPHI(I,LMAX,LP1) ENDDO DO K=2,LP1 DO I=1,IMAX EMPL(I,K)=QH2O(I,K-1)*P(I,K)*(P(I,K)-PRESS(I,K-1))*GP0INV ENDDO ENDDO DO K=LP2,LL DO I=1,IMAX EMPL(I,K)=QH2O(I,K-LMAX)*P(I,K-LMAX)* & (PRESS(I,K-LMAX)-P(I,K-LMAX))*GP0INV ENDDO ENDDO DO I=1,IMAX EMPL(I,LLP1)=EMPL(I,LL) ENDDO DO I=1,IMAX TPL(I,1)=TEMP(I,LMAX) TPL(I,LP1)=0.50*(TTTT(I,LP1)+TEMP(I,LMAX)) TPL(I,LLP1)=0.50*(TTTT(I,LMAX)+TEMP(I,LMAX)) ENDDO DO K=2,LMAX DO I=1,IMAX TPL(I,K)=TTTT(I,K) ENDDO ENDDO DO K=2,LMAX DO I=1,IMAX TPL(I,K+LMAX)=TTTT(I,K) ENDDO ENDDO !---THE CALCULATION OF TTEMP HAS BEEN MOVED TO E1E288 TO SAVE STORAGE DO K=2,LMAX DO I=1,IMAX AVEPHI(I,K,K)=EMXX(I,K-1) ENDDO ENDDO DO I=1,IMAX AVEPHI(I,LP1,LMAX)=AVEPHI(I,LP1,LMAX)+EMPL(I,LMAX) ENDDO ! COMPUTE LOGS OF WATER MASS ARGUMENTS FOR E1E288; ! THE CORRESPONDING QUANTITY FOR E3V88 (EMPL) MUST BE OBTAINED ! WITHIN THAT SUBROUTINE, AS EMPL IS USED AFTER E3V88 IS CALLED. DO KK=1,LP1 DO K=1,LP1 DO I=1,IMAX AVEPHI(I,K,KK)=LOG10(AVEPHI(I,K,KK)) ENDDO ENDDO ENDDO ! !----------------------------------------------------------------------- ! ! CALL E1E288 FOR EMISSIVITY TRANSMISSION FCTNS FOR H2O ! ----------------------------------------------------- CALL E1E288 (E1CTS1,E1CTS2,E1FLX,E1CTW1,E1CTW2,FXO,DT,TEMP) ! !----------------------------------------------------------------------- ! ! CALL E3V88 FOR NBL H2O TRANSMISSIVITIES ! --------------------------------------- CALL E3V88 (EMD,TPL,EMPL) ! !----------------------------------------------------------------------- ! ! COMPUTE NEARBY LAYER AND SPECIAL-CASE TRANSMISSIVITIES FOR EMISS ! USING METHODS FOR H2O GIVEN IN REF. (4) DO K=2,LMAX DO I=1,IMAX CEVAL(I,K-1)=EMD(I,K+LMAX)+EMD(I,K) ENDDO ENDDO DO K=1,LM1 DO I=1,IMAX EMISS(I,LP1,K)=0.50*(EMISS(I,K+1,K+1)+EMISS(I,LP1,K)) ENDDO ENDDO DO K=2,LMAX DO I=1,IMAX EMISS(I,K,K)=CEVAL(I,K-1) ENDDO ENDDO DO I=1,IMAX EMISS(I,LMAX,LP1)=(EMD(I,1)*EMPL(I,1)-EMD(I,LP1)*EMPL(I,LP1))/ & EMX1(I) + 0.25*(EMISS(I,LMAX,LP1)+EMISS(I,LP1,LMAX)) EMISS(I,LP1,LP1)=2.0*EMD(I,LP1) EMISS(I,LP1,LMAX)=2.0*(EMD(I,1)*EMPL(I,1)-EMD(I,LLP1)* & EMPL(I,LLP1))/(QH2O(I,LMAX)*PRESS(I,LMAX)* & (P(I,LP1)-PRESS(I,LMAX))*GP0INV) ENDDO !----------------------------------------------------------------------- ! SUBROUTINE SPA88 IS CALLED TO OBTAIN EXACT CTS FOR WATER, ! CO2 AND O3, AND APPROXIMATE CTS CO2 AND O3 CALCULATIONS. ! CALL SPA88 (PRESS,TEMP) ! !----------------------------------------------------------------------- ! ! THIS SECTION PERFORMS THE CALCULATION OF "EMISSIVITY" ! FLUXES FOR THE 4 COMPONENTS COMPRISING THE LW FREQUENCY REGION. ! (EMISS = THE 0-160,1200-2200 CM-1 BAND; EMISS2 THE 800-990, ! 1070-1200 CM-1 BAND; TO3 THE 990-1070 CM-1 BAND; CO21 THE 560-800 ! CM-1 BAND). EMISST IS THE COMBINED EXCHANGE TERM AND FLX THE ! COMBINED NET FLUX. ! !----------------------------------------------------------------------- DO K=1,LP1 DO I=1,IMAX TC(I,K)=TEMP(I,K)*TEMP(I,K)*TEMP(I,K)*TEMP(I,K) ENDDO ENDDO DO K=2,LP1 DO I=1,IMAX OSS(I,K)=OSOUR(I,K)-OSOUR(I,K-1) CSS(I,K)=CSOUR(I,K)-CSOUR(I,K-1) DTC(I,K)=TC(I,K)-TC(I,K-1) SS2(I,K)=SS1(I,K)-SS1(I,K-1) ENDDO ENDDO !1 DO K=1,LP1 !1 DO I=1,IMAX !1 EMISST(I,1,K)=TC(I,1)*E1FLX(I,K)+SS1(I,1)*EMISS2(I,K,1)+ !1 & OSOUR(I,1)*TO3(I,K,1)+CSOUR(I,1)*CO2SP(I,K) !1 ENDDO !1 ENDDO !********************************************************** !1 DO K=1,LP1 !1 DO KP=2,LP1 !1 DO I=1,IMAX !1 EMISST(I,KP,K)=DTC(I,KP)*EMISS(I,KP,K)+ !1 & SS2(I,KP)*EMISS2(I,KP,K)+OSS(I,KP)*TO3(I,KP,K)+ !1 & CSS(I,KP)*CO21(I,KP,K) !1 ENDDO !1 ENDDO !1 ENDDO !1 DO K=1,LP1 !1 DO KP=1,LP1 !1 DO I=1,IMAX !1 EMISST(I,KP,K)=EMISST(I,KP,K)*CLDFAC(I,KP,K) !1 ENDDO !1 ENDDO !1 ENDDO !---THE ORDER OF THE 911 LOOPS IS TO MINIMIZE "MEMORY CONFLICTS" ! ON THE CYBER 205 DO K=1,LP1 DO I=1,IMAX !f FLX(I,K)=0.0 FLX(I,K)= & (TC(I,1)*E1FLX(I,K)+SS1(I,1)*EMISS2(I,K,1)+ & OSOUR(I,1)*TO3(I,K,1)+CSOUR(I,1)*CO2SP(I,K)) * & cldfac(i,1,k) ENDDO ENDDO !f DO KP=1,LP1 DO KP=2,LP1 DO K=1,LP1 DO I=1,IMAX !f FLX(I,K)=FLX(I,K)+EMISST(I,KP,K) FLX(I,K)=FLX(I,K)+ & ( DTC(I,KP)*EMISS(I,KP,K)+ & SS2(I,KP)*EMISS2(I,KP,K)+OSS(I,KP)*TO3(I,KP,K)+ & CSS(I,KP)*CO21(I,KP,K) ) *cldfac(i,kp,k) ENDDO ENDDO ENDDO !----------------------------------------------------------------------- ! THIS SECTION COMPUTES THE EMISSIVITY CTS HEATING RATES FOR 2 ! EMISSIVITY BANDS: THE 0-160,1200-2200 CM-1 BAND AND THE 800- ! 990,1070-1200 CM-1 BAND. THE REMAINING CTS COMTRIBUTIONS ARE ! CONTAINED IN CTSO3, COMPUTED IN SPA88. !----------------------------------------------------------------------- DO K=1,LMAX DO I=1,IMAX CTS(I,K)=RADCON*DELP(I,K)*(TC(I,K)* & (E1CTW2(I,K)*CLDFAC(I,K+1,1)-E1CTW1(I,K)*CLDFAC(I,K,1)) + & SS1(I,K)*(EMISS2(I,K+1,1)*CLDFAC(I,K+1,1)- & EMISS2(I,K,1)*CLDFAC(I,K,1))) ENDDO ENDDO DO K=1,LMAX DO I=1,IMAX CEVAL(I,K)=TC(I,K)*(CLDFAC(I,K,1)*(E1CTS1(I,K)-E1CTW1(I,K)) - & CLDFAC(I,K+1,1)*(E1CTS2(I,K)-E1CTW2(I,K))) ENDDO ENDDO DO I=1,IMAX FLX1E1(I)=0.0 ENDDO DO K=1,LMAX DO I=1,IMAX FLX1E1(I)=FLX1E1(I)+CEVAL(I,K) ENDDO ENDDO DO I=1,IMAX FLX1E1(I)=FLX1E1(I)+TC(I,LP1)*CLDFAC(I,LP1,1)* & (E1CTS1(I,LP1)-E1CTW1(I,LP1)) ENDDO !----------------------------------------------------------------------- ! FINAL SECTION OBTAINS EMISSIVITY HEATING RATES, ! TOTAL HEATING RATES AND THE FLUX AT THE GROUND !----------------------------------------------------------------------- ! .....CALCULATE THE EMISSIVITY HEATING RATES DO K=1,LMAX DO I=1,IMAX HEATEM(I,K)=RADCON*(FLX(I,K+1)-FLX(I,K))*DELP(I,K) ENDDO ENDDO ! .....CALCULATE THE TOTAL HEATING RATES DO K=1,LMAX DO I=1,IMAX HEATRA(I,K)=HEATEM(I,K)-CTS(I,K)-CTSO3(I,K)+EXCTS(I,K) ENDDO ENDDO ! .....CALCULATE THE FLUX AT EACH FLUX LEVEL USING THE FLUX AT THE ! TOP (FLX1E1+GXCTS) AND THE INTEGRAL OF THE HEATING RATES (VSUM1) DO K=1,LMAX DO I=1,IMAX VSUM1(I,K)=HEATRA(I,K)*DELP2(I,K)*RADCON1 ENDDO ENDDO DO I=1,IMAX TOPFLX(I)=FLX1E1(I)+GXCTS(I) ENDDO !---ONLY THE SURFACE VALUE OF FLUX (GRNFLX) IS NEEDED UNLESS ! THE THICK CLOUD SECTION IS INVOKED. DO I=1,IMAX FLXNET(I,1)=TOPFLX(I) ENDDO DO K=2,LP1 DO I=1,IMAX FLXNET(I,K)=FLXNET(I,K-1)+VSUM1(I,K-1) ENDDO ENDDO IF (Present(LSFC)) THEN DO I=1,IMAX GRNFLX(I)=FLXNET(I,LSFC(I)+1) ENDDO ELSE DO I=1,IMAX GRNFLX(I)=FLXNET(I,LP1) ENDDO ENDIF !----------------------------------------------------------------------- ! THIS IS THE THICK CLOUD SECTION.OPTIONALLY,IF THICK CLOUD ! FLUXES ARE TO BE "CONVECTIVELY ADJUSTED",IE,DF/DP IS CONSTANT, ! FOR CLOUDY PART OF GRID POINT, THE FOLLOWING CODE IS EXECUTED. !***FIRST,COUNT THE NUMBER OF CLOUDS ALONG THE LAT. ROW. SKIP THE ! ENTIRE THICK CLOUD COMPUTATION OF THERE ARE NO CLOUDS. !----------------------------------------------------------------------- ICNT=0 DO I=1,IMAX ICNT=ICNT+NCLDS(I) ENDDO IF (ICNT == 0) GO TO 6999 !***THE OUTER LOOP IS ON THE NO. OF CLOUDS( BEGINNING AT NC=2,AS NC=1 IS ! A DUMMY LAYER) ! ************************************** DO 6001 NC=2,LP1 ! ************************************** !***OBTAIN THE LARGEST NO. OF CLOUDS IN THE LAT. ROW. IF THE CLOUD ! INDEX (NC) IS GREATER THAT THIS NO.,TERMINATE THE CLOUD LOOP ICNT=0 DO I=1,IMAX IF (NCLDS(I) >= NC-1) THEN ICNT=ICNT+1 ENDIF ENDDO IF (ICNT == 0) GO TO 6999 !***FOR THE NC''TH CLOUD LEVEL,OBTAIN THE TOP AND BOTTOM INDEX,AND ! COUNT THE NO. OF THICK CLOUDS (DEFINED AS A NONZERO DIFFERENCE ! BET. THE TOP AND BOTTOM INDICES). SKIP TO THE NEXT CLOUD INDEX ! (INCREMENT NC) IF THIS NO. IS ZERO. DO I=1,IMAX ITOP(I)=KTOP(I,NC) IBOT(I)=KBTM(I,NC) ENDDO ICNT=0 DO I=1,IMAX IF (ITOP(I) < IBOT(I)) THEN ICNT=ICNT+1 ENDIF ENDDO IF (ICNT == 0) GO TO 6001 !***OBTAIN THE PRESSURES AND FLUXES OF THE TOP AND BOTTOM OF ! THE NC'TH CLOUD (IT IS ASSUMED THAT ALL KTOP AND KBTM'S HAVE ! BEEN DEFINED!). DO I=1,IMAX PTOP(I)=P(I,ITOP(I)) FTOP(I)=FLXNET(I,ITOP(I)) PBOT(I)=P(I,IBOT(I)+1) FBOT(I)=FLXNET(I,IBOT(I)+1) ENDDO !***OBTAIN THE "FLUX DERIVATIVE" DF/DP (DELPTC) DO I=1,IMAX DELPTC(I)=(FBOT(I)-FTOP(I))/(PBOT(I)-PTOP(I)) ENDDO !***CALCULATE THE TOT. FLUX CHG. FROM THE TOP OF THE CLOUD, FOR ! ALL LEVELS. THIS QUANTITY WILL BE USED ONLY IN THE THICK CLOUD ! LEVELS, BUT, FOR EFFICIENCY, IS COMPUTED FOR ALL LEVELS. DO K=1,LP1 DO I=1,IMAX Z1(I,K)=(P(I,K)-PTOP(I))*DELPTC(I)+FTOP(I) ENDDO ENDDO !***USING THIS FLUX CHG. IN THE CLOUDY PART OF THE GRID BOX, OBTAIN ! THE NEW FLUXES, WEIGHTING THE CLEAR AND CLOUDY FLUXES:AGAIN, ONLY ! THE FLUXES IN THICK-CLOUD LEVELS WILL EVENTUALLY BE USED. DO K=1,LP1 DO I=1,IMAX FLXTHK(I,K)=FLXNET(I,K)*(1.0-CAMT(I,NC)*EMCLD(I,NC)) + & Z1(I,K)*CAMT(I,NC)*EMCLD(I,NC) ENDDO ENDDO !***MERGE FLXTHK INTO FLXNET FOR APPROPRIATE LEVELS. DO K=1,LP1 DO I=1,IMAX IF (K > ITOP(I) .and. K <= IBOT(I) & .and. (NC-1) <= NCLDS(I)) THEN FLXNET(I,K)=FLXTHK(I,K) ENDIF ENDDO ENDDO ! ******END OF CLOUD LOOP***** 6001 CONTINUE ! *************************** !----------------------------------------------------------------------- ! THE THICK CLOUD SECTION ENDS HERE. !----------------------------------------------------------------------- 6999 CONTINUE !----------------------------------------------------------------------- !***THE FINAL STEP IS TO RECOMPUTE THE HEATING RATES BASED ON THE ! REVISED FLUXES: DO K=1,LMAX DO I=1,IMAX HEATRA(I,K)=RADCON*(FLXNET(I,K+1)-FLXNET(I,K))*DELP(I,K) ENDDO ENDDO !---------------------------------------------------------------------** END SUBROUTINE FST88 !####################################################################### !####################################################################### SUBROUTINE E1E288 (G1,G2,G3,G4,G5,FXOE1,DTE1,TEMP) !----------------------------------------------------------------------- ! SUBROUTINE E1E288 COMPUTES THE EXCHANGE TERMS IN THE FLUX EQUATION ! FOR LONGWAVE RADIATION FOR ALL TERMS EXCEPT THE EXCHANGE WITH THE ! TOP OF THE ATMOSPHERE. THE METHOD IS A TABLE LOOKUP ON A PRE- ! COMPUTED E2 FUNCTION (DEFINED IN REF. (4)). ! THE E1 FUNCTION CALCULATIONS (FORMERLY DONE IN SUBROUTINE ! E1V88 COMPUTE THE FLUX RESULTING FROM THE EXCHANGE OF PHOTONS ! BETWEEN A LAYER AND THE TOP OF THE ATMOSPHERE. THE METHOD IS A ! TABLE LOOKUP ON A PRE-COMPUTED E1 FUNCTION. ! CALCULATIONS ARE DONE IN TWO FREQUENCY RANGES: ! 1) 0-560,1200-2200 CM-1 FOR Q(APPROX) ! 2) 160-560 CM-1 FOR Q(APPROX,CTS). ! MOTIVATION FOR THESE CALCULATIONS IS IN REFERENCES (1) AND (4). !----------------------------------------------------------------------- ! INPUTS: (COMMON BLOCKS) ! TABLE1,TABLE2,TABLE,EM1,EM1WDE TABCOM ! AVEPHI TFCOM ! TEMP ARGUMENT LIST ! T KDACOM ! FXOE1,DTE1 ARGUMENT LIST ! ! OUTPUTS: ! EMISS TFCOM ! G1,G2,G3 ARGUMENT LIST,FOR 1ST FREQ. RANGE ! G4,G5 ARGUMENT LIST,FOR 2ND FREQ. RANGE ! ! CALLED BY : FST88 ! CALLS : ! !----------------------------------------------------------------------- !#include "hcon.h" !----------------------------------------------------------------------- !#include "tabcom.h" real :: EM1 (28,180),EM1WDE(28,180),EM3 (28,180), & & TABLE1(28,180),TABLE2(28,180),TABLE3(28,180) COMMON /TABCOM/ EM1 , EM1WDE, EM3 , & & TABLE1, TABLE2, TABLE3 !----------------------------------------------------------------------- REAL, INTENT(OUT), DIMENSION(:,:) :: G1,G2,G3,G4,G5 REAL, INTENT(IN) , DIMENSION(:,:) :: FXOE1,DTE1,TEMP !**** I/O ARGUMENTS !CC DIMENSION FXOE1(IMAX,LP1),DTE1(IMAX,LP1),G1(IMAX,LP1), !CC & G2(IMAX,LMAX),G3(IMAX,LP1),G4(IMAX,LP1),G5(IMAX,LMAX) !CC DIMENSION TEMP(IMAX,LP1) !----------------------------------------------------------------------- integer :: kk, k, i, kp integer :: IT1(IMAX,LL3P) real :: TTEMP(IMAX,LP1,LP1) real :: FIT1(IMAX,LL3P),DTOT1(IMAX,LL3P), & DTOT2(IMAX,LL3P),DTOT3(IMAX,LL3P),DTOT4(IMAX,LL3P), & WW1(IMAX,LP1),WW2(IMAX,LP1) real :: FYO(IMAX,LP1),DU(IMAX,LP1),F1(IMAX,LP1) !----------------------------------------------------------------------- real :: FXO(IMAX,LP1),FIVAL(IMAX,LP1), & DT(IMAX,LP1),F2(IMAX,LP1),F3(IMAX,LP1) integer :: IVAL(IMAX,LP1) real :: T1(5040),T2(5040),T4(5040) real :: EM1V(5040),EM1VW(5040) ! DIMENSION TTMP1D(IMAX,LP1M) CHARACTER(LEN=40) :: err_string !----------------------------------------------------------------------- EQUIVALENCE (EM1V(1),EM1(1,1)),(EM1VW(1),EM1WDE(1,1)) ! EQUIVALENCE (TTMP1D,TTEMP) ! EQUIVALENCE (IVAL,IT1),(FIVAL,FIT1),(DT,WW1) ! EQUIVALENCE (FXO,WW2),(F2,DTOT2),(F3,DTOT3) EQUIVALENCE (T1(1),TABLE1(1,1)),(T2(1),TABLE2(1,1)), & (T4(1),TABLE3(1,1)) !----------------------------------------------------------------------- DO KK=1,LP1 DO K=1,LP1 DO I=1,IMAX TTEMP(I,K,KK)=TTTT(I,K) ENDDO ENDDO ENDDO DO K=1,LM1 DO I=1,IMAX TTEMP(I,K+1,K+1)=TEMP(I,LMAX) ENDDO ENDDO DO I=1,IMAX TTEMP(I,LP1,LMAX)=TEMP(I,LM1) ENDDO !----------------------------------------------------------------------- !---THE REST OF THE SUBROUTINE IS CARRIED OUT WITHIN A K-LOOP (121) ! PRIMARILY TO CONSERVE STORAGE. SOME LOSS OF SPEED OCCURS. !----------------------------------------------------------------------- DO 121 K=1,LP1 !----------------------------------------------------------------------- ! DO KP=1,LP1 DO I=1,IMAX FXO(I,KP)=AINT((TTEMP(I,KP,K)-100.)*0.10+1.0) FYO(I,KP)=AINT((AVEPHI(I,KP,K)+16.)*10.) ENDDO ENDDO DO KP=1,LP1 DO I=1,IMAX DU(I,KP)=AVEPHI(I,KP,K)-0.10*(FYO(I,KP)+1.0)+16.1 DT(I,KP)=TTEMP(I,KP,K)-10.*FXO(I,KP)-90. ENDDO ENDDO DO KP=1,LP1 DO I=1,IMAX FYO(I,KP)=28.*FYO(I,KP) ENDDO ENDDO DO KP=1,LP1 DO I=1,IMAX FIVAL(I,KP)=FYO(I,KP)+FXO(I,KP) ENDDO ENDDO DO KP=1,LP1 DO I=1,IMAX IVAL(I,KP)=FIVAL(I,KP) ENDDO ENDDO DO KP=1,LP1 DO I=1,IMAX IF (IVAL(I,KP) > 5040) THEN Write (err_string(1:32),800) i,kp,ival(i,kp) 800 Format ('i,kp,ival=',2i6,i10) CALL Error_Mesg ('E1E288 in LONGWAVE_MOD', & err_string(1:32), NOTE) ENDIF F1(I,KP)=T1(IVAL(I,KP)) F2(I,KP)=T2(IVAL(I,KP)) F3(I,KP)=T4(IVAL(I,KP)) ENDDO ENDDO DO KP=1,LP1 DO I=1,IMAX EMISS(I,KP,K)=F1(I,KP)+DU(I,KP)*F2(I,KP)+DT(I,KP)*F3(I,KP) ENDDO ENDDO !***THE FOLLOWING IS THE CALCULATION FOR THE E1 FUNCTION, FORMERLY ! DONE IN SUBROUTINE E1V88. THE MOVE TO E1E288 IS DUE TO THE ! SAVINGS IN OBTAINING INDEX VALUES (THE TEMP. INDICES HAVE ! BEEN OBTAINED IN FST88, WHILE THE U-INDICES ARE OBTAINED ! IN THE E2 CALCS.,WITH K=1). ! ! ---------------- IF (K == 1) THEN ! ---------------- DO KP=1,LP1 DO I=1,IMAX FIT1(I,KP)=FYO(I,KP)+FXOE1(I,KP) ENDDO ENDDO DO KP=1,LMAX DO I=1,IMAX FIT1(I,LP1+KP)=FYO(I,KP+1)+FXOE1(I,KP) ENDDO ENDDO DO KP=1,LP1 DO I=1,IMAX FIT1(I,KP+LLP1)=FYO(I,KP)+FXOE1(I,1) ENDDO ENDDO DO KP=1,LL3P DO I=1,IMAX IT1(I,KP)=FIT1(I,KP) ENDDO ENDDO DO KP=1,LP1 DO I=1,IMAX WW1(I,KP)=10.-DTE1(I,KP) ENDDO ENDDO DO KP=1,LP1 DO I=1,IMAX WW2(I,KP)=0.10-DU(I,KP) ENDDO ENDDO DO KP=1,LL3P DO I=1,IMAX DTOT1(I,KP)=EM1V(IT1(I,KP)) DTOT2(I,KP)=EM1V(IT1(I,KP)+1) DTOT3(I,KP)=EM1V(IT1(I,KP)+28) DTOT4(I,KP)=EM1V(IT1(I,KP)+29) ENDDO ENDDO DO KP=1,LP1 DO I=1,IMAX G1(I,KP)=WW1(I,KP)*WW2(I,KP)*DTOT1(I,KP)+ & WW2(I,KP)*DTE1(I,KP)*DTOT2(I,KP)+ & WW1(I,KP)*DU(I,KP)*DTOT3(I,KP)+ & DTE1(I,KP)*DU(I,KP)*DTOT4(I,KP) ENDDO ENDDO DO KP=1,LMAX DO I=1,IMAX G2(I,KP)=WW1(I,KP)*WW2(I,KP+1)*DTOT1(I,KP+LP1)+ & WW2(I,KP+1)*DTE1(I,KP)*DTOT2(I,KP+LP1)+ & WW1(I,KP)*DU(I,KP+1)*DTOT3(I,KP+LP1)+ & DTE1(I,KP)*DU(I,KP+1)*DTOT4(I,KP+LP1) ENDDO ENDDO DO KP=1,LP1 DO I=1,IMAX G3(I,KP)=WW1(I,1)*WW2(I,KP)*DTOT1(I,KP+LLP1)+ & WW2(I,KP)*DTE1(I,1)*DTOT2(I,KP+LLP1)+ & WW1(I,1)*DU(I,KP)*DTOT3(I,KP+LLP1)+ & DTE1(I,1)*DU(I,KP)*DTOT4(I,KP+LLP1) ENDDO ENDDO !...REPEAT FOR WIDE-BAND COMPUTATIONS OF G4,G5 DO KP=1,LLP1 DO I=1,IMAX DTOT1(I,KP)=EM1VW(IT1(I,KP)) DTOT2(I,KP)=EM1VW(IT1(I,KP)+1) DTOT3(I,KP)=EM1VW(IT1(I,KP)+28) DTOT4(I,KP)=EM1VW(IT1(I,KP)+29) ENDDO ENDDO DO KP=1,LP1 DO I=1,IMAX G4(I,KP)=WW1(I,KP)*WW2(I,KP)*DTOT1(I,KP)+ & WW2(I,KP)*DTE1(I,KP)*DTOT2(I,KP)+ & WW1(I,KP)*DU(I,KP)*DTOT3(I,KP)+ & DTE1(I,KP)*DU(I,KP)*DTOT4(I,KP) ENDDO ENDDO DO KP=1,LMAX DO I=1,IMAX G5(I,KP)=WW1(I,KP)*WW2(I,KP+1)*DTOT1(I,KP+LP1)+ & WW2(I,KP+1)*DTE1(I,KP)*DTOT2(I,KP+LP1)+ & WW1(I,KP)*DU(I,KP+1)*DTOT3(I,KP+LP1)+ & DTE1(I,KP)*DU(I,KP+1)*DTOT4(I,KP+LP1) ENDDO ENDDO ! ------------- ENDIF !----------------------------------------------------------------------- 121 CONTINUE !----------------------------------------------------------------------- END SUBROUTINE E1E288 !####################################################################### !####################################################################### SUBROUTINE E3V88 (EMV,TV,AV) !C IMPLICIT NONE !----------------------------------------------------------------------- ! SUBROUTINE E3V88 COMPUTES NEARBY LAYER TRANSMISSIVITIES FOR ! H2O USING A TABLE LOOKUP OF THE PRE-COMPUTED E3 FUNCTION ! ( DESCRIBED IN REF. (4)). !----------------------------------------------------------------------- ! INPUTS: (COMMON BLOCKS,ARGS.) ! TV,AV ARGUMENT LIST ! EM3 TABCOM ! ! OUTPUTS: ! EMV ARGUMENT LIST ! ! CALLED BY : FST88 ! !----------------------------------------------------------------------- !#include "hcon.h" !----------------------------------------------------------------------- !#include "tabcom.h" real :: EM1 (28,180),EM1WDE(28,180),EM3 (28,180), & TABLE1(28,180),TABLE2(28,180),TABLE3(28,180) COMMON /TABCOM/ EM1 ,EM1WDE ,EM3 , & TABLE1 ,TABLE2 ,TABLE3 !----------------------------------------------------------------------- REAL, INTENT(OUT), DIMENSION(:,:) :: EMV REAL, INTENT(IN), DIMENSION(:,:) :: TV,AV INTEGER, DIMENSION(SIZE(TV,1),SIZE(TV,2)) :: IT REAL, DIMENSION(SIZE(TV,1),SIZE(TV,2)) :: FXO,FYO,FIT,WW1, & AV1LOG,D1,D2,D3,D4,TVAL,DT,WW2,UVAL,DU !----------------------------------------------------------------------- REAL :: EM3V(5040) INTEGER I,K !----------------------------------------------------------------------- EQUIVALENCE (EM3V(1),EM3(1,1)) !----------------------------------------------------------------------- AV1LOG(:,:)=LOG10(AV(:,:)) FXO(:,:)=AINT((TV(:,:)-100.)*0.10+1.0) TVAL(:,:)=90.+10.*FXO(:,:) DT(:,:)=TV(:,:)-TVAL(:,:) FYO(:,:)=AINT((AV1LOG(:,:)+16.)*10.) UVAL(:,:)=-16.1+0.10*(FYO(:,:)+1.0) DU(:,:)=AV1LOG(:,:)-UVAL(:,:) FIT(:,:)=FXO(:,:)+FYO(:,:)*28. IT(:,:)=FIT(:,:) DO K=1,SIZE(TV,2) DO I=1,SIZE(TV,1) D1(I,K)=EM3V(IT(I,K)) D2(I,K)=EM3V(IT(I,K)+1) D3(I,K)=EM3V(IT(I,K)+28) D4(I,K)=EM3V(IT(I,K)+29) ENDDO ENDDO WW1(:,:)=10.-DT(:,:) WW2(:,:)=0.10-DU(:,:) EMV(:,:)=WW1(:,:)*WW2(:,:)*D1(:,:)+WW2(:,:)*DT(:,:)*D2(:,:)+ & WW1(:,:)*DU(:,:)*D3(:,:)+DT(:,:)*DU(:,:)*D4(:,:) !----------------------------------------------------------------------- END SUBROUTINE E3V88 !####################################################################### !####################################################################### SUBROUTINE SPA88 (PRESS,TEMP) !----------------------------------------------------------------------- ! SUBROUTINE SPA88 COMPUTES EXACT CTS HEATING RATES AND FLUXES AND ! CORRESPONDING CTS EMISSIVITY QUANTITIES FOR H2O,CO2 AND O3. ! ! INPUTS: (COMMON BLOCKS) ! ACOMB,BCOMB,APCM,BPCM BDCOMB ! ATPCM,BTPCM,BETACM BDCOMB ! BETINW BDWIDE ! VAR1,VAR2,P,DELP,DELP2 KDACOM ! TOTVO2,TO3SPC,CO2SP1,CO2SP2,CO2SP TFCOM ! CLDFAC CLDCOM ! SKO2D TABCOM ! SORC,CSOUR,OSOUR SRCCOM ! OUTPUTS: ! EXCTS,EXCTSN,CTSO3 TFCOM ! GXCTS,FCTSG RDFLUX ! CALLED BY: ! FST88 ! !----------------------------------------------------------------------- !#include "hcon.h" !----------------------------------------------------------------------- !#include "bandta.h" !----------------------------------------------------------------------- !#include "bdwide.h" !----------------------------------------------------------------------- !#include "bdcomb.h" !----------------------------------------------------------------------- REAL, INTENT(IN), DIMENSION(:,:) :: PRESS,TEMP ! DIMENSION PRESS(IMAX,LP1),TEMP(IMAX,LP1) !----------------------------------------------------------------------- integer :: k, n real :: CAPPHI(IMAX,LP1),CAPPSI(IMAX,LP1),TT(IMAX,LMAX), & FAC1(IMAX,LMAX),FAC2(IMAX,LMAX), & CTMP(IMAX,LP1),CDIF(IMAX,LMAX), & TOPM(IMAX,LMAX),TOPPHI(IMAX,LMAX), & BOT(IMAX), & X(IMAX,LMAX),Y(IMAX,LMAX),VFLUX(IMAX,LMAX) !----------------------------------------------------------------------- real :: F(IMAX,LMAX),FF(IMAX,LMAX),AG(IMAX,LMAX),AGG(IMAX,LMAX), & AH(IMAX,LMAX),AHH(IMAX,LMAX),AI(IMAX,LMAX),AII(IMAX,LMAX), & AJ(IMAX,LMAX),AJJ(IMAX,LMAX),FELSX(IMAX,LMAX), & PHITMP(IMAX,LMAX),PSITMP(IMAX,LMAX) !----------------------------------------------------------------------- ! EQUIVALENCE (CAPPHI,AJ,AI,AH,AG,F,PHITMP) ! EQUIVALENCE (CAPPSI,AJJ,AII,AHH,AGG,FF,PSITMP) !----------------------------------------------------------------------- ! COMPUTE TEMP. TEMPEDENCE OF LINE INTENSITIES (CAPPHI,CAPPSI) DO K=1,LMAX DO I=1,IMAX X(I,K)=TEMP(I,K)-250. ENDDO ENDDO DO K=1,LMAX DO I=1,IMAX Y(I,K)=X(I,K)*X(I,K) ENDDO ENDDO DO K=1,LMAX DO I=1,IMAX EXCTS(I,K)=0.0 ENDDO ENDDO DO I=1,IMAX GXCTS(I)=0.0 ENDDO !*************BEGIN LOOP ON FREQUENCY BANDS (N)************************* DO 200 N=1,NBLM !*********************************************************************** DO K=1,LMAX DO I=1,IMAX F(I,K)=0.044194*(APCM(N)*X(I,K)+BPCM(N)*Y(I,K)) FF(I,K)=0.044194*(ATPCM(N)*X(I,K)+BTPCM(N)*Y(I,K)) ENDDO ENDDO DO K=1,LMAX DO I=1,IMAX AG(I,K)=(1.418191+F(I,K))*F(I,K)+1.0 AGG(I,K)=(1.418191+FF(I,K))*FF(I,K)+1.0 ENDDO ENDDO DO K=1,LMAX DO I=1,IMAX AH(I,K)=AG(I,K)*AG(I,K) AHH(I,K)=AGG(I,K)*AGG(I,K) ENDDO ENDDO DO K=1,LMAX DO I=1,IMAX AI(I,K)=AH(I,K)*AH(I,K) AII(I,K)=AHH(I,K)*AHH(I,K) ENDDO ENDDO DO K=1,LMAX DO I=1,IMAX AJ(I,K)=AI(I,K)*AI(I,K) AJJ(I,K)=AII(I,K)*AII(I,K) ENDDO ENDDO DO K=1,LMAX DO I=1,IMAX CAPPHI(I,K)=AJ(I,K)*AJ(I,K) CAPPSI(I,K)=AJJ(I,K)*AJJ(I,K) ENDDO ENDDO ! OBTAIN WEIGHED OPTICAL PATH AND MEAN PRESSURE (TOPM,TOPPHI) DO K=1,LMAX DO I=1,IMAX PHITMP(I,K)=VAR1(I,K)*CAPPHI(I,K) PSITMP(I,K)=VAR2(I,K)*CAPPSI(I,K) ENDDO ENDDO DO I=1,IMAX TOPM(I,1)=PHITMP(I,1) TOPPHI(I,1)=PSITMP(I,1) ENDDO DO K=2,LMAX ! DO I=1,IMAX TOPM(:,K)=TOPM(:,K-1)+PHITMP(:,K) TOPPHI(:,K)=TOPPHI(:,K-1)+PSITMP(:,K) ! ENDDO ENDDO !***COMPUTE CTS TRANSMISSION FCTNS FOR ALL BANDS (TT) DO K=1,LMAX DO I=1,IMAX FAC1(I,K)=ACOMB(N)*TOPM(I,K) ENDDO ENDDO DO K=1,LMAX DO I=1,IMAX FAC2(I,K)=FAC1(I,K)*TOPM(I,K)/(BCOMB(N)*TOPPHI(I,K)) ENDDO ENDDO DO K=1,LMAX DO I=1,IMAX FELSX(I,K)=FAC1(I,K)/SQRT(1.0+FAC2(I,K)) ENDDO ENDDO !***ADD H2O CONTINUUM PATHS FOR FREQ. BANDS 5-14 IF (N >= 5 .and. N <= 14) THEN DO K=1,LMAX DO I=1,IMAX FELSX(I,K)=FELSX(I,K)+BETACM(N)*TOTVO2(I,K+1)*SKO2D ENDDO ENDDO ENDIF !***ADD O3 PATHS FOR BAND 13:990-1070 BAND IF (N == 13) THEN DO K=1,LMAX DO I=1,IMAX FELSX(I,K)=FELSX(I,K)+TO3SPC(I,K) ENDDO ENDDO ENDIF DO K=1,LMAX DO I=1,IMAX FELSX(I,K)=-1.0*FELSX(I,K) ENDDO ENDDO !***COMPUTE EXPONENTIAL (TT) DO K=1,LMAX DO I=1,IMAX TT(I,K)=EXP(FELSX(I,K)) ENDDO ENDDO !***ADD THE CO2 TRANSMISSIVITIES FOR BANDS 9 AND 10 (15 UM BAND) IF (N == 9) THEN DO K=1,LMAX DO I=1,IMAX TT(I,K)=TT(I,K)*CO2SP1(I,K+1) ENDDO ENDDO ENDIF IF (N == 10) THEN DO K=1,LMAX DO I=1,IMAX TT(I,K)=TT(I,K)*CO2SP2(I,K+1) ENDDO ENDDO ENDIF !*** COMPUTE EXACT CTS HEATING RATES (EXCTS,EXCTSN) DO I=1,IMAX CTMP(I,1)=1.0 ENDDO DO K=2,LP1 DO I=1,IMAX CTMP(I,K)=TT(I,K-1)*CLDFAC(I,K,1) ENDDO ENDDO DO K=1,LMAX DO I=1,IMAX CDIF(I,K)=SORC(I,K,N)*(CTMP(I,K+1)-CTMP(I,K)) ENDDO ENDDO DO K=1,LMAX DO I=1,IMAX EXCTSN(I,K,N)=RADCON*DELP(I,K)*CDIF(I,K) ENDDO ENDDO DO K=1,LMAX DO I=1,IMAX EXCTS(I,K)=EXCTS(I,K)+EXCTSN(I,K,N) ENDDO ENDDO ! OBTAIN EXACT CTS FLUXES (GXCTS,FCTSG) AND APPROX CTS HEATING ! RATES (CTSO3) DO I=1,IMAX BOT(I)=0.50*DELP(I,LMAX)*(TT(I,LM1)*(P(I,LP1)-PRESS(I,LMAX)) + & TT(I,LMAX)*(P(I,LP1)+PRESS(I,LMAX)-2.0*P(I,LMAX))) ENDDO DO I=1,IMAX FCTSG(I,N)=CLDFAC(I,LP1,1)*(TT(I,LMAX)*SORC(I,LMAX,N)+ & BOT(I)*(SORC(I,LP1,N)-SORC(I,LMAX,N))) ENDDO DO I=1,IMAX GXCTS(I)=GXCTS(I)+FCTSG(I,N) ENDDO !**************END LOOP ON FREQUENCY BANDS****************************** 200 CONTINUE !*********************************************************************** DO K=1,LMAX DO I=1,IMAX VFLUX(I,K)=EXCTS(I,K)*DELP2(I,K)*RADCON1 ENDDO ENDDO DO K=1,LMAX DO I=1,IMAX GXCTS(I)=GXCTS(I)-VFLUX(I,K) ENDDO ENDDO !*** COMPUTE APPROXIMATE CTS HEATING RATES FOR 15UM AND 9.6 UM BANDS ! (CTSO3) DO K=1,LMAX DO I=1,IMAX CTSO3(I,K)=RADCON*DELP(I,K)* & (CSOUR(I,K)*(CO2SP(I,K+1)*CLDFAC(I,K+1,1)- & CO2SP(I,K)*CLDFAC(I,K,1)) + & OSOUR(I,K)*(TO3(I,K+1,1)*CLDFAC(I,K+1,1)- & TO3(I,K,1)*CLDFAC(I,K,1))) ENDDO ENDDO !----------------------------------------------------------------------- END SUBROUTINE SPA88 !####################################################################### !####################################################################### ! CO21 TFCOM SUBROUTINE NLTE (PRESS,RRVCO2) !----------------------------------------------------------------------- ! SUBROUTINE NLTE IS THE PRESENT FORMULATION OF AN NLTE ! CALCULATION OF THE SOURCE FCTN IN THE 15 UM REGION (2 BANDS) !----------------------------------------------------------------------- ! ! RRVCO2 IS THE VOLUME MIXING RATIO OF CO2 (SCALAR) ! !----------------------------------------------------------------------- ! ! INPUTS (COMMON BLOCKS) ! BANDLO,BANDHI BANDTA ! DELP KDACOM ! CO21 TFCOM ! CSOUR1,CSOUR2 SRCCOM ! OUTPUTS: ! SORC SRCCOM !----------------------------------------------------------------------- REAL, INTENT(IN), DIMENSION(:,:) :: PRESS REAL, INTENT(IN) :: RRVCO2 !----------------------------------------------------------------------- !#include "hcon.h" !----------------------------------------------------------------------- !#include "bandta.h" !----------------------------------------------------------------------- integer :: n, k, i, KP, IITER real :: CMTRX(IMAX,LMAX,INLTE),VNL(IMAX,LMAX,INLTE), & VSUM(IMAX,INLTE),BDENOM(IMAX,INLTE),CDIAG(IMAX,INLTE), & TCOLL(IMAX,INLTE),FNLTE(IMAX,INLTE), & AZ(IMAX,INLTE),AG(IMAX,INLTE), & C1B7(2),C2B7(2),CENT(2),DEL(2) !----------------------------------------------------------------------- !*****DEGENERACY FACTOR=0.5***** real DEGEN DATA DEGEN /0.5/ !----------------------------------------------------------------------- DO N=1,2 CENT(N)=0.50*(BANDLO(N+NNLTE)+BANDHI(N+NNLTE)) DEL(N)=BANDHI(N+NNLTE)-BANDLO(N+NNLTE) C1B7(N)=3.7412E-5*CENT(N)*CENT(N)*CENT(N)*DEL(N) C2B7(N)=1.4387*CENT(N) ENDDO DO K=1,INLTE DO I=1,IMAX TCOLL(I,K)=DEGEN*1.5E-5*PRESS(I,LP1)/(SECPDA*PRESS(I,K)) FNLTE(I,K)=3.5*TCOLL(I,K)*C1B7(1)/(RRVCO2*C2B7(1)) ENDDO ENDDO DO K=1,INLTE DO KP=1,LM1 DO I=1,IMAX CMTRX(I,KP,K)=RADCON*DELP(I,K)*(CO21(I,KP,K+1)-CO21(I,KP+1,K+1) & -CO21(I,KP,K)+CO21(I,KP+1,K)) ENDDO ENDDO ENDDO DO K=1,INLTE DO I=1,IMAX CMTRX(I,LMAX,K)=DELP(I,K)*(CO21(I,LMAX,K+1)-CO21(I,LMAX,K))* & RADCON ENDDO ENDDO DO K=1,INLTE DO I=1,IMAX CDIAG(I,K)=CMTRX(I,K,K) ENDDO ENDDO DO K=1,INLTE DO I=1,IMAX BDENOM(I,K)=1.0/(1.0-FNLTE(I,K)*CDIAG(I,K)) SORC(I,K,9)=CSOUR1(I,K) ENDDO ENDDO DO K=1,INLTE DO KP=1,LMAX DO I=1,IMAX VNL(I,KP,K)=CMTRX(I,KP,K)*CSOUR1(I,KP) ENDDO ENDDO ENDDO DO K=1,INLTE DO I=1,IMAX VSUM(I,K)=0.0 ENDDO ENDDO DO K=1,INLTE DO KP=INLTEP,LMAX DO I=1,IMAX VSUM(I,K)=VSUM(I,K)+VNL(I,KP,K) ENDDO ENDDO ENDDO DO K=1,INLTE DO I=1,IMAX AZ(I,K)=CSOUR1(I,K)+FNLTE(I,K)*VSUM(I,K) ENDDO ENDDO ! **** ITERATION BEGINS HERE **** DO IITER=1,2 ! ------------ DO K=1,INLTE DO KP=1,INLTE DO I=1,IMAX VNL(I,KP,K)=CMTRX(I,KP,K)*SORC(I,KP,9) ENDDO ENDDO ENDDO DO K=1,INLTE DO I=1,IMAX VSUM(I,K)=0.0 ENDDO ENDDO DO K=1,INLTE DO KP=1,INLTE DO I=1,IMAX VSUM(I,K)=VSUM(I,K)+VNL(I,KP,K) ENDDO ENDDO ENDDO DO K=1,INLTE DO I=1,IMAX AG(I,K)=FNLTE(I,K)*(VSUM(I,K)-CDIAG(I,K)*SORC(I,K,9)) SORC(I,K,9)=BDENOM(I,K)*(AZ(I,K)+AG(I,K)) ENDDO ENDDO ! ----------- ENDDO ! ----------- DO K=1,INLTE DO I=1,IMAX TCOLL(I,K)=DEGEN*1.5E-5*PRESS(I,LP1)/(SECPDA*PRESS(I,K)) FNLTE(I,K)=3.5*TCOLL(I,K)*C1B7(2)/(RRVCO2*C2B7(2)) ENDDO ENDDO DO K=1,INLTE DO KP=1,LM1 DO I=1,IMAX CMTRX(I,KP,K)=RADCON*DELP(I,K)*(CO21(I,KP,K+1)-CO21(I,KP+1,K+1) & -CO21(I,KP,K)+CO21(I,KP+1,K)) ENDDO ENDDO ENDDO DO K=1,INLTE DO I=1,IMAX CMTRX(I,LMAX,K)=DELP(I,K)*(CO21(I,LMAX,K+1)-CO21(I,LMAX,K))* & RADCON ENDDO ENDDO DO K=1,INLTE DO I=1,IMAX CDIAG(I,K)=CMTRX(I,K,K) ENDDO ENDDO DO K=1,INLTE DO I=1,IMAX BDENOM(I,K)=1.0/(1.0-FNLTE(I,K)*CDIAG(I,K)) SORC(I,K,10)=CSOUR2(I,K) ENDDO ENDDO DO K=1,INLTE DO KP=1,LMAX DO I=1,IMAX VNL(I,KP,K)=CMTRX(I,KP,K)*CSOUR2(I,KP) ENDDO ENDDO ENDDO DO K=1,INLTE DO I=1,IMAX VSUM(I,K)=0.0 ENDDO ENDDO DO K=1,INLTE DO KP=INLTEP,LMAX DO I=1,IMAX VSUM(I,K)=VSUM(I,K)+VNL(I,KP,K) ENDDO ENDDO ENDDO DO K=1,INLTE DO I=1,IMAX AZ(I,K)=CSOUR2(I,K)+FNLTE(I,K)*VSUM(I,K) ENDDO ENDDO ! **** ITERATION BEGINS HERE **** DO IITER=1,2 ! ------------ DO K=1,INLTE DO KP=1,INLTE DO I=1,IMAX VNL(I,KP,K)=CMTRX(I,KP,K)*SORC(I,KP,10) ENDDO ENDDO ENDDO DO K=1,INLTE DO I=1,IMAX VSUM(I,K)=0.0 ENDDO ENDDO DO K=1,INLTE DO KP=1,INLTE DO I=1,IMAX VSUM(I,K)=VSUM(I,K)+VNL(I,KP,K) ENDDO ENDDO ENDDO DO K=1,INLTE DO I=1,IMAX AG(I,K)=FNLTE(I,K)*(VSUM(I,K)-CDIAG(I,K)*SORC(I,K,10)) SORC(I,K,10)=BDENOM(I,K)*(AZ(I,K)+AG(I,K)) ENDDO ENDDO ! ----------- ENDDO ! ----------- END SUBROUTINE NLTE !####################################################################### !####################################################################### SUBROUTINE CLO88 (KTOP,KBTM,NCLDS,CAMT,EMCLD) !----------------------------------------------------------------------- !#include "hcon.h" !----------------------------------------------------------------------- ! SUBROUTINE CLO88 COMPUTES CLOUD TRANSMISSION FUNCTIONS FOR THE ! LONGWAVE CODE,USING CODE WRITTEN BY BERT KATZ (301-763-8161). ! AND MODIFIED BY DAN SCHWARZKOPF IN DECEMBER,1988. !----------------------------------------------------------------------- ! INPUT PARAMETERS ! ---------------- ! ! KTOP = INDEX OF (DATA LEVEL) PRESSURE OF CLOUD TOP,USED ! IN THE LONGWAVE PROGRAM ! KBTM = INDEX OF (DATA LEVEL) PRESSURE OF CLOUD BOTTOM, ! USED IN THE LONGWAVE PROGRAM ! NCLDS = NO. CLOUDS AT EACH GRID PT. ! EMCLD = CLOUD EMISSIVITY. SET TO ONE BY DEFAULT, BUT MAY ! BE MODIFIED FOR USE IN LONGWAVE PROGRAM. ! CAMT = CLOUD AMOUNTS OF CLOUDS (THEIR LOCATIONS ARE ! SPECIFIED IN THE KTOP/KBTM INDICES) ! !----------------------------------------------------------------------- ! ! OUTPUT: ! CLDFAC CLDCOM ! !----------------INPUT ARGUMENTS---------------------------------------- INTEGER, INTENT(IN), DIMENSION(:,:) :: KTOP,KBTM INTEGER, INTENT(IN), DIMENSION(:) :: NCLDS REAL, INTENT(IN), DIMENSION(:,:) :: EMCLD,CAMT ! D I M E N S I O N ! & KTOP(IMAX,LP1),KBTM(IMAX,LP1),NCLDS(IMAX),EMCLD(IMAX,LP1), ! & CAMT(IMAX,LP1) !----------------------------------------------------------------------- ! DIMENSION V1(IMAX,LP1,LP1), integer :: k, kp, nc, icnt real :: CAMTV(IMAX),V11D(IMAX,LP1) !----------------------------------------------------------------------- LOGICAL BCLDUD,BCLDSA,BCLDSB,BCLDSC DIMENSION BCLDUD(2*IMAX,LP1),BCLDSA(IMAX,LP1), & BCLDSB(IMAX,LP1),BCLDSC(IMAX,LP1) !----------------------------------------------------------------------- ! ***LOGICAL ARRAY DEFINITIONS*** ! ! BCLDUD = T IF (AT GRID POINT I) :1) VERTICAL LEVEL K IS BELOW ! THE TOP OF CLOUD (NC);2) VERTICAL LEVEL K IS NOT ! BELOW THE BOTTOM OF CLOUD (NC). FOR EACH VERTICAL ! LEVEL K,THE FIRST (IMAX) PTS. IS FOR TEST (1); ! THE SECOND (IMAX) PTS. IS FOR TEST (2). ! BCLDSA = T IF (AT GRID POINT I) EITHER VERTICAL LEVEL K OR ! VERTICAL LEVEL KP IS BELOW THE TOP OF CLOUD (NC). ! BCLDSB = T IF (AT GRID POINT I) EITHER VERTICAL LEVEL K OR ! VERTICAL LEVEL KP IS NOT BELOW THE BOTTOM OF ! CLOUD (NC). ! BCLDSC = T IF (AT GRID POINT I) :1)VERTICAL LEVEL K IS ABOVE CLOUD ! (NC) AND VERTICAL LEVEL KP IS BELOW CLOUD (NC), OR VICE ! VERSA; 2) VERTICAL LEVELS K,KP ARE BOTH BELOW THE TOP ! OF CLOUD (NC) AND ABOVE THE BOTTOM OF CLOUD (NC) ! BELOW THE TOP OF CLOUD (NC). !----------------------------------------------------------------------- ! *** Initialize CLDFAC *** ! *** Note: If the user wants to initialize this quantity ! *** then it must be passed as an argument (in & out) DO K =1,LP1 DO KP=1,LP1 DO I =1,IMAX CLDFAC(I,KP,K)=1.0 ENDDO ENDDO ENDDO ! ******* LOOP ON CLOUDS ******* DO 1 NC=1,LMAX ! ****************************** ! *** IF NO GRID PT HAS AN (NCTH) CLOUD, *** ! *** THE CLOUD CALCULATION IS COMPLETE *** ICNT=0 DO I=1,IMAX IF (NC <= NCLDS(I)) THEN ICNT=ICNT+1 ENDIF ENDDO IF (ICNT == 0) then RETURN endif ! ****************************** !***COMPUTE BCLDUD*** DO K=1,LP1 DO I=1,IMAX BCLDUD(I,K)=K.GT.KTOP(I,NC+1) BCLDUD(IMAX+I,K)=K.LE.KBTM(I,NC+1) ENDDO ENDDO !***INITIALIZE THE CLOUD TF FOR THE (NCTH) CLOUD TO UNITY. ! DO K=1,LP1 ! DO KP=1,LP1 ! DO I=1,IMAX ! V1(I,KP,K)=1.0 ! ENDDO ! ENDDO ! ENDDO !***DEFINE THE CLOUD TRANSMISSIVITY. THE TRANSMISSIVITY IS TAKEN AS ! (1-CLOUD EMISSIVITY) FOR THE FRACTION OF THE GRID BOX WITH THE NCTH ! LAYER OF CLOUD, AND ONE OTHERWISE. ( NOTE THAT IF AT POINT I, THE ! NCTH CLOUD DOES NOT EXIST, ITS AMOUNT IS TAKEN AS ZERO). AS A RESULT, ! THE CLOUD TRANSMISSIVITY MAY BE DEFINED AS THOUGH THE EMISSIVITY ! WERE UNITY AND THE CLOUD AMOUNT WERE MULTIPLIED BY THE EMISSIVITY. DO I=1,IMAX IF (NC <= NCLDS(I)) THEN CAMTV(I)=1.0-CAMT(I,NC+1)*EMCLD(I,NC+1) ELSE CAMTV(I)=1.0 ENDIF ENDDO DO K=1,LP1 DO I=1,IMAX V11D(I,K)=CAMTV(I) ENDDO ENDDO !***THE 221 LOOP IS THE K-LOOP. WE USED SUCH A LOOP PRIMARILY TO REDUCE ! SPACE*** DO 221 K=1,LP1 DO KP=1,LP1 DO I=1,IMAX BCLDSA(I,KP)=BCLDUD(I,KP).or.BCLDUD(I,K) BCLDSB(I,KP)=BCLDUD(IMAX+I,KP).or.BCLDUD(IMAX+I,K) ENDDO ENDDO !---COMBINE BCLDSA AND BCLDSB INTO BCLDSC DO KP=1,LP1 DO I=1,IMAX BCLDSC(I,KP)=BCLDSA(I,KP).and.BCLDSB(I,KP) ENDDO ENDDO !***DEFINE TRANSMISSION FUNCTION FOR LEVELS WHERE LOGICAL ARRAY BCLDSC ! INDICATES THAT A CLOUD IS PRESENT ! DO KP=1,LP1 ! DO I=1,IMAX ! IF (BCLDSC(I,KP)) THEN ! V1(I,KP,K)=V11D(I,KP) ! ENDIF ! ENDDO ! ENDDO !221 CONTINUE !***RANDOM OVERLAP OF CLOUD IS ASSUMED IS OBTAINING THE TOTAL CLOUD TF** ! DO K=1,LP1 DO KP=1,LP1 DO I=1,IMAX ! CLDFAC(I,KP,K)=CLDFAC(I,KP,K)*V1(I,KP,K) if (bcldsc(i,kp) ) then CLDFAC(I,KP,K)=CLDFAC(I,KP,K)*V11d(I,KP) else endif ENDDO ENDDO ! ENDDO 221 CONTINUE ! ****************************** 1 CONTINUE ! ****************************** !----------------------------------------------------------------------- END SUBROUTINE CLO88 !####################################################################### !####################################################################### SUBROUTINE TABLE !----------------------------------------------------------------------- ! SUBROUTINE TABLE COMPUTES TABLE ENTRIES USED IN THE LONGWAVE RADIA ! PROGRAM. ALSO CALCULATED ARE INDICES USED IN STRIP-MINING AND FOR ! SOME PRE-COMPUTABLE FUNCTIONS. !----------------------------------------------------------------------- !#include "hcon.h" !----------------------------------------------------------------------- !#include "bandta.h" !----------------------------------------------------------------------- !#include "bdwide.h" !----------------------------------------------------------------------- !#include "bdcomb.h" !----------------------------------------------------------------------- !#include "tabcom.h" integer :: i1, i2e, i2, indx, jp, n, icnt, j, ia, nsubds, nsb real :: cent, del, bdlo, bdhi, c1, anu real :: EM1 (28,180),EM1WDE(28,180),EM3 (28,180), & TABLE1(28,180),TABLE2(28,180),TABLE3(28,180) COMMON /TABCOM/ EM1 ,EM1WDE ,EM3 , & TABLE1 ,TABLE2 ,TABLE3 !----------------------------------------------------------------------- real :: SUM(28,180),PERTSM(28,180),SUM3(28,180),SUMWDE(28,180), & SRCWD(28,NBLX),SRC1NB(28,NBLW),DBDTNB(28,NBLW) real :: ZMASS(181),ZROOT(181),SC(28),DSC(28),XTEMV(28), & TFOUR(28),FORTCU(28),X(28),X1(28),X2(180),SRCS(28), & SUM4(28),SUM6(28),SUM7(28),SUM8(28),SUM4WD(28), & R1(28),R2(28),S2(28),T3(28),R1WD(28) real :: EXPO(180),FAC(180) real :: CNUSB(30),DNUSB(30) real :: ALFANB(NBLW),AROTNB(NBLW) real :: ANB(NBLW),BNB(NBLW),CENTNB(NBLW),DELNB(NBLW), & BETANB(NBLW) ! !*** NOTE: THE DATA,EQUIVALENCE AND DIMENSION STATEMENTS FOR QUANTITIES ! EQUIVALENCED TO COMMON BLOCK BANDTA DEPEND ON THE VALUE OF THE ! PARAMETER NBLW. ! !----------------------------------------------------------------------- !***COMPUTE LOCAL QUANTITIES AND AO3,BO3,AB15 !....FOR NARROW-BANDS... DO N=1,NBLW ANB(N)=ARNDM(N) BNB(N)=BRNDM(N) CENTNB(N)=0.50*(BANDLO(N)+BANDHI(N)) DELNB(N)=BANDHI(N)-BANDLO(N) BETANB(N)=BETAD(N) ENDDO AB15(1)=ANB(57)*BNB(57) AB15(2)=ANB(58)*BNB(58) !....FOR WIDE BANDS... AB15WD=AWIDE*BWIDE !***COMPUTE INDICES: INDX2,KMAXV ICNT=0 DO I1=1,LMAX I2E=LP1-I1 DO I2=1,I2E ICNT=ICNT+1 INDX=LP1*(I2-1)+LP2*I1 INDX1(ICNT)=MOD(INDX,LP1) IF (INDX1(ICNT) == 0) INDX1(ICNT)=LP1 INDX2(ICNT)=(INDX+LMAX)/LP1 ENDDO ENDDO KMAXV(1)=1 DO I=2,LMAX KMAXV(I)=KMAXV(I-1)+(LP2-I) ENDDO KMAXVM=KMAXV(LMAX) !***COMPUTE RATIOS OF CONT. COEFFS SKC1R=BETAWD/BETINW SKO3R=BETAD(61)/BETINW SKO2D=1.0/BETINW ! !****BEGIN TABLE COMPUTATIONS HERE*** !***COMPUTE TEMPS, MASSES FOR TABLE ENTRIES !---NOTE: THE DIMENSIONING AND INITIALIZATION OF XTEMV AND OTHER ARRAYS ! WITH DIMENSION OF 28 IMPLY A RESTRICTION OF MODEL TEMPERATURES FROM ! 100K TO 370K. !---THE DIMENSIONING OF ZMASS,ZROOT AND OTHER ARRAYS WITH DIMENSION OF ! 180 IMPLY A RESTRICTION OF MODEL H2O AMOUNTS SUCH THAT OPTICAL PATHS ! ARE BETWEEN 10**-16 AND 10**2, IN CGS UNITS. ZMASS(1)=1.0E-16 DO J=1,180 JP=J+1 ZROOT(J)=SQRT(ZMASS(J)) ZMASS(JP)=ZMASS(J)*1.258925411 ENDDO DO I=1,28 XTEMV(I)=90.+10.*I TFOUR(I)=XTEMV(I)*XTEMV(I)*XTEMV(I)*XTEMV(I) FORTCU(I)=4.0*XTEMV(I)*XTEMV(I)*XTEMV(I) ENDDO !******THE COMPUTATION OF SOURCE,DSRCE IS NEEDED ONLY ! FOR THE COMBINED WIDE-BAND CASE.TO OBTAIN THEM,THE SOURCE ! MUST BE COMPUTED FOR EACH OF THE (NBLX) WIDE BANDS(=SRCWD) ! THEN COMBINED (USING IBAND) INTO SOURCE. DO N=1,NBLY DO I=1,28 SOURCE(I,N)=0.0 ENDDO ENDDO DO N=1,NBLX DO I=1,28 SRCWD(I,N)=0.0 ENDDO ENDDO !---BEGIN FREQ. LOOP (ON N) DO 211 N=1,NBLX !***THE 160-1200 BAND CASES IF (N <= 46) THEN CENT=CENTNB(N+16) DEL=DELNB(N+16) BDLO=BANDLO(N+16) BDHI=BANDHI(N+16) ENDIF !***THE 2270-2380 BAND CASE IF (N == NBLX) THEN CENT=CENTNB(NBLW) DEL=DELNB(NBLW) BDLO=BANDLO(NBLW) BDHI=BANDHI(NBLW) ENDIF !***FOR PURPOSES OF ACCURACY, ALL EVALUATIONS OF PLANCK FCTNS ARE MADE ! ON 10 CM-1 INTERVALS, THEN SUMMED INTO THE (NBLX) WIDE BANDS. NSUBDS=(DEL-1.0E-1)/10+1 DO 213 NSB=1,NSUBDS IF (NSB .ne. NSUBDS) THEN CNUSB(NSB)=10.*(NSB-1)+BDLO+5.0 DNUSB(NSB)=10. ELSE CNUSB(NSB)=0.50*(10.*(NSB-1)+BDLO+BDHI) DNUSB(NSB)=BDHI-(10.*(NSB-1)+BDLO) ENDIF C1=(3.7412E-5)*CNUSB(NSB)**3 !---BEGIN TEMP. LOOP (ON I) DO I=1,28 X(I)=1.4387*CNUSB(NSB)/XTEMV(I) X1(I)=EXP(X(I)) SRCS(I)=C1/(X1(I)-1.0) SRCWD(I,N)=SRCWD(I,N)+SRCS(I)*DNUSB(NSB) ENDDO 213 CONTINUE 211 CONTINUE !***THE FOLLOWING LOOPS CREATE THE COMBINED WIDE BAND QUANTITIES SOURCE ! AND DSRCE DO N=1,40 DO I=1,28 SOURCE(I,IBAND(N))=SOURCE(I,IBAND(N))+SRCWD(I,N) ENDDO ENDDO DO N=9,NBLY DO I=1,28 SOURCE(I,N)=SRCWD(I,N+32) ENDDO ENDDO DO N=1,NBLY DO I=1,27 DSRCE(I,N)=(SOURCE(I+1,N)-SOURCE(I,N))*0.10 ENDDO ENDDO DO N=1,NBLW ALFANB(N)=BNB(N)*ANB(N) AROTNB(N)=SQRT(ALFANB(N)) ENDDO !***FIRST COMPUTE PLANCK FCTNS (SRC1NB) AND DERIVATIVES (DBDTNB) FOR ! USE IN TABLE EVALUATIONS. THESE ARE DIFFERENT FROM SOURCE,DSRCE ! BECAUSE DIFFERENT FREQUENCY PTS ARE USED IN EVALUATION, THE FREQ. ! RANGES ARE DIFFERENT, AND THE DERIVATIVE ALGORITHM IS DIFFERENT. ! DO 301 N=1,NBLW CENT=CENTNB(N) DEL=DELNB(N) !---NOTE: AT PRESENT, THE IA LOOP IS ONLY USED FOR IA=2. THE LOOP STRUCT ! IS KEPT SO THAT IN THE FUTURE, WE MAY USE A QUADRATURE SCHEME FOR ! THE PLANCK FCTN EVALUATION, RATHER THAN USE THE MID-BAND FREQUENCY. DO IA=1,3 ANU=CENT+0.50*(IA-2)*DEL C1=(3.7412E-5)*ANU*ANU*ANU+1.E-20 !---TEMPERATURE LOOP--- DO I=1,28 X(I)=1.4387*ANU/XTEMV(I) X1(I)=EXP(X(I)) SC(I)=C1/((X1(I)-1.0)+1.E-20) DSC(I)=SC(I)*SC(I)*X(I)*X1(I)/(XTEMV(I)*C1) ENDDO IF (IA == 2) THEN DO I=1,28 SRC1NB(I,N)=DEL*SC(I) DBDTNB(I,N)=DEL*DSC(I) ENDDO ENDIF ENDDO 301 CONTINUE !***NEXT COMPUTE R1,R2,S2,AND T3- COEFFICIENTS USED FOR E3 FUNCTION ! WHEN THE OPTICAL PATH IS LESS THAN 10-4. IN THIS CASE, WE ASSUME A ! DIFFERENT DEPENDENCE ON (ZMASS). !---ALSO OBTAIN R1WD, WHICH IS R1 SUMMED OVER THE 160-560 CM-1 RANGE DO I=1,28 SUM4(I)=0.0 SUM6(I)=0.0 SUM7(I)=0.0 SUM8(I)=0.0 SUM4WD(I)=0.0 ENDDO DO N=1,NBLW CENT=CENTNB(N) !***PERFORM SUMMATIONS FOR FREQ. RANGES OF 0-560,1200-2200 CM-1 FOR SUM4 ! SUM6,SUM7,SUM8 IF (CENT < 560. .or. CENT > 1200. .and. CENT <= 2200.) THEN DO I=1,28 SUM4(I)=SUM4(I)+SRC1NB(I,N) SUM6(I)=SUM6(I)+DBDTNB(I,N) SUM7(I)=SUM7(I)+DBDTNB(I,N)*AROTNB(N) SUM8(I)=SUM8(I)+DBDTNB(I,N)*ALFANB(N) ENDDO ENDIF !***PERFORM SUMMATIONS OVER 160-560 CM-1 FREQ RANGE FOR E1 CALCS (SUM4WD IF (CENT > 160. .and. CENT < 560.) THEN DO I=1,28 SUM4WD(I)=SUM4WD(I)+SRC1NB(I,N) ENDDO ENDIF ENDDO DO I=1,28 R1(I)=SUM4(I)/TFOUR(I) R2(I)=SUM6(I)/FORTCU(I) S2(I)=SUM7(I)/FORTCU(I) T3(I)=SUM8(I)/FORTCU(I) R1WD(I)=SUM4WD(I)/TFOUR(I) ENDDO DO J=1,180 DO I=1,28 SUM(I,J)=0.0 PERTSM(I,J)=0.0 SUM3(I,J)=0.0 SUMWDE(I,J)=0.0 ENDDO ENDDO !------- FREQUENCY LOOP BEGINS ----------------------------------------- DO 411 N=1,NBLW CENT=CENTNB(N) !***PERFORM CALCULATIONS FOR FREQ. RANGES OF 0-560,1200-2200 CM-1 IF (CENT < 560. .or. CENT > 1200. .and. CENT <= 2200.) THEN DO J=1,180 X2(J)=AROTNB(N)*ZROOT(J) EXPO(J)=EXP(-X2(J)) ENDDO DO J=1,180 IF (X2(J) >= 100.) THEN EXPO(J)=0.0 ENDIF ENDDO DO J=121,180 FAC(J)=ZMASS(J)*(1.0-(1.0+X2(J))*EXPO(J))/(X2(J)*X2(J)) ENDDO DO J=1,180 DO I=1,28 SUM(I,J)=SUM(I,J)+SRC1NB(I,N)*EXPO(J) PERTSM(I,J)=PERTSM(I,J)+DBDTNB(I,N)*EXPO(J) ENDDO ENDDO DO J=121,180 DO I=1,28 SUM3(I,J)=SUM3(I,J)+DBDTNB(I,N)*FAC(J) ENDDO ENDDO ENDIF !---COMPUTE SUM OVER 160-560 CM-1 RANGE FOR USE IN E1 CALCS (SUMWDE) IF (CENT > 160. .and. CENT < 560.) THEN DO J=1,180 DO I=1,28 SUMWDE(I,J)=SUMWDE(I,J)+SRC1NB(I,N)*EXPO(J) ENDDO ENDDO ENDIF 411 CONTINUE !----------------------------------------------------------------------- DO J=1,180 DO I=1,28 EM1(I,J)=SUM(I,J)/TFOUR(I) TABLE1(I,J)=PERTSM(I,J)/FORTCU(I) ENDDO ENDDO DO J=121,180 DO I=1,28 EM3(I,J)=SUM3(I,J)/FORTCU(I) ENDDO ENDDO DO J=1,179 DO I=1,28 TABLE2(I,J)=(TABLE1(I,J+1)-TABLE1(I,J))*10. ENDDO ENDDO DO J=1,180 DO I=1,27 TABLE3(I,J)=(TABLE1(I+1,J)-TABLE1(I,J))*0.10 ENDDO ENDDO DO I=1,28 TABLE2(I,180)=0.0 ENDDO DO J=1,180 TABLE3(28,J)=0.0 ENDDO DO J=1,2 DO I=1,28 EM1(I,J)=R1(I) ENDDO ENDDO DO J=1,120 DO I=1,28 EM3(I,J)=R2(I)/2.0-S2(I)*SQRT(ZMASS(J))/3.0+T3(I)*ZMASS(J)/8.0 ENDDO ENDDO DO J=121,180 DO I=1,28 EM3(I,J)=EM3(I,J)/ZMASS(J) ENDDO ENDDO !***NOW COMPUTE E1 TABLES FOR 160-560 CM-1 BANDS ONLY. ! WE USE R1WD AND SUMWDE OBTAINED ABOVE. DO J=1,180 DO I=1,28 EM1WDE(I,J)=SUMWDE(I,J)/TFOUR(I) ENDDO ENDDO DO J=1,2 DO I=1,28 EM1WDE(I,J)=R1WD(I) ENDDO ENDDO !----------------------------------------------------------------------- END SUBROUTINE TABLE !####################################################################### SUBROUTINE RAD_ALLOC (IMAX) INTEGER,INTENT(IN) :: IMAX !----------------------------------------------------------------------- ! Allocates temporary space for longwave code !----------------------------------------------------------------------- !-------------- SRCCOM ----------------- ALLOCATE (SORC(IMAX,LP1,NBLY)) ALLOCATE (CSOUR1(IMAX,LP1)) ALLOCATE (CSOUR2(IMAX,LP1)) ALLOCATE (OSOUR(IMAX,LP1)) ALLOCATE (CSOUR(IMAX,LP1)) ALLOCATE (SS1(IMAX,LP1)) !-------------- KDACOM ----------------- ALLOCATE (QH2O(IMAX,LP1)) ALLOCATE (P(IMAX,LP1)) ALLOCATE (DELP2(IMAX,LMAX)) ALLOCATE (DELP(IMAX,LMAX)) ALLOCATE (TTTT(IMAX,LP1)) ALLOCATE (VAR1(IMAX,LMAX)) ALLOCATE (VAR2(IMAX,LMAX)) ALLOCATE (VAR3(IMAX,LMAX)) ALLOCATE (VAR4(IMAX,LMAX)) ALLOCATE (CNTVAL(IMAX,LP1)) !-------------- RDFLUX ----------------- ALLOCATE (FLX1E1(IMAX)) ALLOCATE (GXCTS(IMAX)) ALLOCATE (FCTSG(IMAX,NBLY)) !-------------- CLDCOM ----------------- ALLOCATE (CLDFAC(IMAX,LP1,LP1)) !-------------- TFCOM ----------------- ALLOCATE (TO3 (IMAX,LP1,LP1)) ALLOCATE (CO21 (IMAX,LP1,LP1)) ALLOCATE (EMISS (IMAX,LP1,LP1)) ALLOCATE (EMISS2(IMAX,LP1,LP1)) ALLOCATE (AVEPHI(IMAX,LP1,LP1)) ALLOCATE (CTS (IMAX,LMAX)) ALLOCATE (CTSO3(IMAX,LMAX)) ALLOCATE (EXCTS(IMAX,LMAX)) ALLOCATE (EXCTSN(IMAX,LMAX,NBLY)) ALLOCATE (E1FLX (IMAX,LP1)) ALLOCATE (CO2NBL(IMAX,LMAX)) ALLOCATE (CO2SP1(IMAX,LP1)) ALLOCATE (CO2SP2(IMAX,LP1)) ALLOCATE (CO2SP (IMAX,LP1)) ALLOCATE (TO3SPC(IMAX,LMAX)) ALLOCATE (TOTVO2(IMAX,LP1)) !----------------------------------------------------------------------- !-------------- TABCOM + Initialize tables for longwave ---------------- IF (IMAX .ne. IMAXold .or. LMAX .ne. LMAXold) THEN IF (Allocated(INDX1)) DeAllocate (INDX1) IF (Allocated(INDX2)) DeAllocate (INDX2) IF (Allocated(KMAXV)) DeAllocate (KMAXV) Allocate (INDX1(LP1V)) Allocate (INDX2(LP1V)) Allocate (KMAXV(LP1)) CALL TABLE IMAXold=IMAX LMAXold=LMAX ENDIF !----------------------------------------------------------------------- END SUBROUTINE RAD_ALLOC !####################################################################### SUBROUTINE RAD_DEALLOC !----------------------------------------------------------------------- ! ----- SRCCOM ----- DEALLOCATE (SORC) DEALLOCATE (CSOUR1) DEALLOCATE (CSOUR2) DEALLOCATE (OSOUR) DEALLOCATE (CSOUR) DEALLOCATE (SS1) ! ----- KDACOM ----- DEALLOCATE (QH2O) DEALLOCATE (P) DEALLOCATE (DELP2) DEALLOCATE (DELP) DEALLOCATE (TTTT) DEALLOCATE (VAR1) DEALLOCATE (VAR2) DEALLOCATE (VAR3) DEALLOCATE (VAR4) DEALLOCATE (CNTVAL) ! ----- RDFLUX ----- DEALLOCATE (FLX1E1) DEALLOCATE (GXCTS) DEALLOCATE (FCTSG) ! ----- CLDCOM ----- DEALLOCATE (CLDFAC) !-------------- TABCOM ----------------- ! DEALLOCATE (INDX1) ! DEALLOCATE (INDX2) ! DEALLOCATE (KMAXV) !-------------- TFCOM ----------------- DEALLOCATE (TO3 ) DEALLOCATE (CO21 ) DEALLOCATE (EMISS ) DEALLOCATE (EMISS2) DEALLOCATE (AVEPHI) DEALLOCATE (CTS ) DEALLOCATE (CTSO3) DEALLOCATE (EXCTS) DEALLOCATE (EXCTSN) DEALLOCATE (E1FLX ) DEALLOCATE (CO2NBL) DEALLOCATE (CO2SP1) DEALLOCATE (CO2SP2) DEALLOCATE (CO2SP ) DEALLOCATE (TO3SPC) DEALLOCATE (TOTVO2) !----------------------------------------------------------------------- END SUBROUTINE RAD_DEALLOC !####################################################################### END MODULE LONGWAVE_MOD