| 1 | ! (c) British Crown Copyright 2008, the Met Office. |
|---|
| 2 | ! All rights reserved. |
|---|
| 3 | ! |
|---|
| 4 | ! Redistribution and use in source and binary forms, with or without modification, are permitted |
|---|
| 5 | ! provided that the following conditions are met: |
|---|
| 6 | ! |
|---|
| 7 | ! * Redistributions of source code must retain the above copyright notice, this list |
|---|
| 8 | ! of conditions and the following disclaimer. |
|---|
| 9 | ! * Redistributions in binary form must reproduce the above copyright notice, this list |
|---|
| 10 | ! of conditions and the following disclaimer in the documentation and/or other materials |
|---|
| 11 | ! provided with the distribution. |
|---|
| 12 | ! * Neither the name of the Met Office nor the names of its contributors may be used |
|---|
| 13 | ! to endorse or promote products derived from this software without specific prior written |
|---|
| 14 | ! permission. |
|---|
| 15 | ! |
|---|
| 16 | ! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY EXPRESS OR |
|---|
| 17 | ! IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND |
|---|
| 18 | ! FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR |
|---|
| 19 | ! CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL |
|---|
| 20 | ! DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, |
|---|
| 21 | ! DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER |
|---|
| 22 | ! IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT |
|---|
| 23 | ! OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. |
|---|
| 24 | |
|---|
| 25 | MODULE MOD_COSP |
|---|
| 26 | USE MOD_COSP_TYPES |
|---|
| 27 | USE MOD_COSP_SIMULATOR |
|---|
| 28 | IMPLICIT NONE |
|---|
| 29 | |
|---|
| 30 | CONTAINS |
|---|
| 31 | |
|---|
| 32 | |
|---|
| 33 | !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
|---|
| 34 | !--------------------- SUBROUTINE COSP --------------------------- |
|---|
| 35 | !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
|---|
| 36 | SUBROUTINE COSP(overlap,Ncolumns,cfg,vgrid,gbx,sgx,sgradar,sglidar,isccp,misr,stradar,stlidar) |
|---|
| 37 | |
|---|
| 38 | ! Arguments |
|---|
| 39 | integer,intent(in) :: overlap ! overlap type in SCOPS: 1=max, 2=rand, 3=max/rand |
|---|
| 40 | integer,intent(in) :: Ncolumns |
|---|
| 41 | type(cosp_config),intent(in) :: cfg ! Configuration options |
|---|
| 42 | type(cosp_vgrid),intent(in) :: vgrid ! Information on vertical grid of stats |
|---|
| 43 | type(cosp_gridbox),intent(inout) :: gbx |
|---|
| 44 | type(cosp_subgrid),intent(inout) :: sgx ! Subgrid info |
|---|
| 45 | type(cosp_sgradar),intent(inout) :: sgradar ! Output from radar simulator |
|---|
| 46 | type(cosp_sglidar),intent(inout) :: sglidar ! Output from lidar simulator |
|---|
| 47 | type(cosp_isccp),intent(inout) :: isccp ! Output from ISCCP simulator |
|---|
| 48 | type(cosp_misr),intent(inout) :: misr ! Output from MISR simulator |
|---|
| 49 | type(cosp_radarstats),intent(inout) :: stradar ! Summary statistics from radar simulator |
|---|
| 50 | type(cosp_lidarstats),intent(inout) :: stlidar ! Summary statistics from lidar simulator |
|---|
| 51 | |
|---|
| 52 | ! Local variables |
|---|
| 53 | integer :: Npoints ! Number of gridpoints |
|---|
| 54 | integer :: Nlevels ! Number of levels |
|---|
| 55 | integer :: Nhydro ! Number of hydrometeors |
|---|
| 56 | integer :: Niter ! Number of calls to cosp_simulator |
|---|
| 57 | integer :: i_first,i_last ! First and last gridbox to be processed in each iteration |
|---|
| 58 | integer :: i,j,k,Ni |
|---|
| 59 | integer,dimension(2) :: ix,iy |
|---|
| 60 | logical :: reff_zero |
|---|
| 61 | real :: minv,maxv |
|---|
| 62 | real :: maxp,minp |
|---|
| 63 | integer,dimension(:),allocatable :: & ! Dimensions nPoints |
|---|
| 64 | seed ! It is recommended that the seed is set to a different value for each model |
|---|
| 65 | ! gridbox it is called on, as it is possible that the choice of the same |
|---|
| 66 | ! seed value every time may introduce some statistical bias in the results, |
|---|
| 67 | ! particularly for low values of NCOL. |
|---|
| 68 | ! Types used in one iteration |
|---|
| 69 | type(cosp_gridbox) :: gbx_it |
|---|
| 70 | type(cosp_subgrid) :: sgx_it |
|---|
| 71 | type(cosp_vgrid) :: vgrid_it |
|---|
| 72 | type(cosp_sgradar) :: sgradar_it |
|---|
| 73 | type(cosp_sglidar) :: sglidar_it |
|---|
| 74 | type(cosp_isccp) :: isccp_it |
|---|
| 75 | type(cosp_misr) :: misr_it |
|---|
| 76 | type(cosp_radarstats) :: stradar_it |
|---|
| 77 | type(cosp_lidarstats) :: stlidar_it |
|---|
| 78 | |
|---|
| 79 | !++++++++++ Dimensions ++++++++++++ |
|---|
| 80 | Npoints = gbx%Npoints |
|---|
| 81 | Nlevels = gbx%Nlevels |
|---|
| 82 | Nhydro = gbx%Nhydro |
|---|
| 83 | |
|---|
| 84 | !++++++++++ Apply sanity checks to inputs ++++++++++ |
|---|
| 85 | ! call cosp_check_input('longitude',gbx%longitude,min_val=0.0,max_val=360.0) |
|---|
| 86 | call cosp_check_input('longitude',gbx%longitude,min_val=-180.0,max_val=180.0) |
|---|
| 87 | call cosp_check_input('latitude',gbx%latitude,min_val=-90.0,max_val=90.0) |
|---|
| 88 | call cosp_check_input('dlev',gbx%dlev,min_val=0.0) |
|---|
| 89 | call cosp_check_input('p',gbx%p,min_val=0.0) |
|---|
| 90 | call cosp_check_input('ph',gbx%ph,min_val=0.0) |
|---|
| 91 | call cosp_check_input('T',gbx%T,min_val=0.0) |
|---|
| 92 | call cosp_check_input('q',gbx%q,min_val=0.0) |
|---|
| 93 | call cosp_check_input('sh',gbx%sh,min_val=0.0) |
|---|
| 94 | call cosp_check_input('dtau_s',gbx%dtau_s,min_val=0.0) |
|---|
| 95 | call cosp_check_input('dtau_c',gbx%dtau_c,min_val=0.0) |
|---|
| 96 | call cosp_check_input('dem_s',gbx%dem_s,min_val=0.0,max_val=1.0) |
|---|
| 97 | call cosp_check_input('dem_c',gbx%dem_c,min_val=0.0,max_val=1.0) |
|---|
| 98 | ! Point information (Npoints) |
|---|
| 99 | call cosp_check_input('land',gbx%land,min_val=0.0,max_val=1.0) |
|---|
| 100 | call cosp_check_input('psfc',gbx%psfc,min_val=0.0) |
|---|
| 101 | call cosp_check_input('sunlit',gbx%sunlit,min_val=0.0,max_val=1.0) |
|---|
| 102 | call cosp_check_input('skt',gbx%skt,min_val=0.0) |
|---|
| 103 | ! TOTAL and CONV cloud fraction for SCOPS |
|---|
| 104 | call cosp_check_input('tca',gbx%tca,min_val=0.0,max_val=1.0) |
|---|
| 105 | call cosp_check_input('cca',gbx%cca,min_val=0.0,max_val=1.0) |
|---|
| 106 | ! Precipitation fluxes on model levels |
|---|
| 107 | call cosp_check_input('rain_ls',gbx%rain_ls,min_val=0.0) |
|---|
| 108 | call cosp_check_input('rain_cv',gbx%rain_cv,min_val=0.0) |
|---|
| 109 | call cosp_check_input('snow_ls',gbx%snow_ls,min_val=0.0) |
|---|
| 110 | call cosp_check_input('snow_cv',gbx%snow_cv,min_val=0.0) |
|---|
| 111 | call cosp_check_input('grpl_ls',gbx%grpl_ls,min_val=0.0) |
|---|
| 112 | ! Hydrometeors concentration and distribution parameters |
|---|
| 113 | call cosp_check_input('mr_hydro',gbx%mr_hydro,min_val=0.0) |
|---|
| 114 | ! Effective radius [m]. (Npoints,Nlevels,Nhydro) |
|---|
| 115 | call cosp_check_input('Reff',gbx%Reff,min_val=0.0) |
|---|
| 116 | reff_zero=.true. |
|---|
| 117 | if (any(gbx%Reff > 1.e-8)) then |
|---|
| 118 | reff_zero=.false. |
|---|
| 119 | ! reff_zero == .false. |
|---|
| 120 | ! and gbx%use_reff == .true. Reff use in radar and lidar |
|---|
| 121 | ! and reff_zero == .false. Reff use in lidar and set to 0 for radar |
|---|
| 122 | endif |
|---|
| 123 | if ((gbx%use_reff) .and. (reff_zero)) then ! Inconsistent choice. Want to use Reff but not inputs passed |
|---|
| 124 | print *, '---------- COSP ERROR ------------' |
|---|
| 125 | print *, '' |
|---|
| 126 | print *, 'use_reff==.true. but Reff is always zero' |
|---|
| 127 | print *, '' |
|---|
| 128 | print *, '----------------------------------' |
|---|
| 129 | stop |
|---|
| 130 | endif |
|---|
| 131 | if ((.not. gbx%use_reff) .and. (reff_zero)) then ! No Reff in radar. Default in lidar |
|---|
| 132 | gbx%Reff = DEFAULT_LIDAR_REFF |
|---|
| 133 | print *, '---------- COSP WARNING ------------' |
|---|
| 134 | print *, '' |
|---|
| 135 | print *, 'Using default Reff in lidar simulations' |
|---|
| 136 | print *, '' |
|---|
| 137 | print *, '----------------------------------' |
|---|
| 138 | endif |
|---|
| 139 | |
|---|
| 140 | ! Aerosols concentration and distribution parameters |
|---|
| 141 | call cosp_check_input('conc_aero',gbx%conc_aero,min_val=0.0) |
|---|
| 142 | ! Checks for CRM mode |
|---|
| 143 | if (Ncolumns == 1) then |
|---|
| 144 | if (gbx%use_precipitation_fluxes) then |
|---|
| 145 | print *, '---------- COSP ERROR ------------' |
|---|
| 146 | print *, '' |
|---|
| 147 | print *, 'Use of precipitation fluxes not supported in CRM mode (Ncolumns=1)' |
|---|
| 148 | print *, '' |
|---|
| 149 | print *, '----------------------------------' |
|---|
| 150 | stop |
|---|
| 151 | endif |
|---|
| 152 | if ((maxval(gbx%dtau_c) > 0.0).or.(maxval(gbx%dem_c) > 0.0)) then |
|---|
| 153 | print *, '---------- COSP ERROR ------------' |
|---|
| 154 | print *, '' |
|---|
| 155 | print *, ' dtau_c > 0.0 or dem_c > 0.0. In CRM mode (Ncolumns=1), ' |
|---|
| 156 | print *, ' the optical depth (emmisivity) of all clouds must be ' |
|---|
| 157 | print *, ' passed through dtau_s (dem_s)' |
|---|
| 158 | print *, '' |
|---|
| 159 | print *, '----------------------------------' |
|---|
| 160 | stop |
|---|
| 161 | endif |
|---|
| 162 | endif |
|---|
| 163 | |
|---|
| 164 | |
|---|
| 165 | ! We base the seed in the decimal part of the surface pressure. |
|---|
| 166 | allocate(seed(Npoints)) |
|---|
| 167 | seed = int(gbx%psfc) ! This is to avoid division by zero when Npoints = 1 |
|---|
| 168 | ! Roj Oct/2008 ... Note: seed value of 0 caused me some problems + I want to |
|---|
| 169 | ! randomize for each call to COSP even when Npoints ==1 |
|---|
| 170 | minp = minval(gbx%psfc) |
|---|
| 171 | maxp = maxval(gbx%psfc) |
|---|
| 172 | if (Npoints .gt. 1) seed=int((gbx%psfc-minp)/(maxp-minp)*100000) + 1 |
|---|
| 173 | |
|---|
| 174 | |
|---|
| 175 | if (gbx%Npoints_it >= gbx%Npoints) then ! One iteration gbx%Npoints |
|---|
| 176 | call cosp_iter(overlap,seed,cfg,vgrid,gbx,sgx,sgradar,sglidar,isccp,misr,stradar,stlidar) |
|---|
| 177 | else ! Several iterations to save memory |
|---|
| 178 | Niter = gbx%Npoints/gbx%Npoints_it ! Integer division |
|---|
| 179 | if (Niter*gbx%Npoints_it < gbx%Npoints) Niter = Niter + 1 |
|---|
| 180 | do i=1,Niter |
|---|
| 181 | i_first = (i-1)*gbx%Npoints_it + 1 |
|---|
| 182 | i_last = i_first + gbx%Npoints_it - 1 |
|---|
| 183 | i_last = min(i_last,gbx%Npoints) |
|---|
| 184 | Ni = i_last - i_first + 1 |
|---|
| 185 | if (i == 1) then |
|---|
| 186 | ! Allocate types for all but last iteration |
|---|
| 187 | call construct_cosp_gridbox(gbx%time,gbx%radar_freq,gbx%surface_radar,gbx%use_mie_tables,gbx%use_gas_abs, & |
|---|
| 188 | gbx%do_ray,gbx%melt_lay,gbx%k2,Ni,Nlevels,Ncolumns,N_HYDRO,gbx%Nprmts_max_hydro, & |
|---|
| 189 | gbx%Naero,gbx%Nprmts_max_aero,Ni,gbx%lidar_ice_type,gbx%isccp_top_height, & |
|---|
| 190 | gbx%isccp_top_height_direction,gbx%isccp_overlap,gbx%isccp_emsfc_lw, & |
|---|
| 191 | gbx%use_precipitation_fluxes,gbx%use_reff, & |
|---|
| 192 | gbx%plat,gbx%sat,gbx%inst,gbx%nchan,gbx%ZenAng, & |
|---|
| 193 | gbx%Ichan(1:gbx%nchan),gbx%surfem(1:gbx%nchan), & |
|---|
| 194 | gbx%co2,gbx%ch4,gbx%n2o,gbx%co, & |
|---|
| 195 | gbx_it) |
|---|
| 196 | call construct_cosp_vgrid(gbx_it,vgrid%Nlvgrid,vgrid%use_vgrid,vgrid%csat_vgrid,vgrid_it) |
|---|
| 197 | call construct_cosp_subgrid(Ni, Ncolumns, Nlevels, sgx_it) |
|---|
| 198 | call construct_cosp_sgradar(cfg,Ni,Ncolumns,Nlevels,N_HYDRO,sgradar_it) |
|---|
| 199 | call construct_cosp_sglidar(cfg,Ni,Ncolumns,Nlevels,N_HYDRO,PARASOL_NREFL,sglidar_it) |
|---|
| 200 | call construct_cosp_isccp(cfg,Ni,Ncolumns,Nlevels,isccp_it) |
|---|
| 201 | call construct_cosp_misr(cfg,Ni,misr_it) |
|---|
| 202 | call construct_cosp_radarstats(cfg,Ni,Ncolumns,vgrid%Nlvgrid,N_HYDRO,stradar_it) |
|---|
| 203 | call construct_cosp_lidarstats(cfg,Ni,Ncolumns,vgrid%Nlvgrid,N_HYDRO,PARASOL_NREFL,stlidar_it) |
|---|
| 204 | elseif (i == Niter) then ! last iteration |
|---|
| 205 | call free_cosp_gridbox(gbx_it,.true.) |
|---|
| 206 | call free_cosp_subgrid(sgx_it) |
|---|
| 207 | call free_cosp_vgrid(vgrid_it) |
|---|
| 208 | call free_cosp_sgradar(sgradar_it) |
|---|
| 209 | call free_cosp_sglidar(sglidar_it) |
|---|
| 210 | call free_cosp_isccp(isccp_it) |
|---|
| 211 | call free_cosp_misr(misr_it) |
|---|
| 212 | call free_cosp_radarstats(stradar_it) |
|---|
| 213 | call free_cosp_lidarstats(stlidar_it) |
|---|
| 214 | ! Allocate types for iterations |
|---|
| 215 | call construct_cosp_gridbox(gbx%time,gbx%radar_freq,gbx%surface_radar,gbx%use_mie_tables,gbx%use_gas_abs, & |
|---|
| 216 | gbx%do_ray,gbx%melt_lay,gbx%k2,Ni,Nlevels,Ncolumns,N_HYDRO,gbx%Nprmts_max_hydro, & |
|---|
| 217 | gbx%Naero,gbx%Nprmts_max_aero,Ni,gbx%lidar_ice_type,gbx%isccp_top_height, & |
|---|
| 218 | gbx%isccp_top_height_direction,gbx%isccp_overlap,gbx%isccp_emsfc_lw, & |
|---|
| 219 | gbx%use_precipitation_fluxes,gbx%use_reff, & |
|---|
| 220 | gbx%plat,gbx%sat,gbx%inst,gbx%nchan,gbx%ZenAng, & |
|---|
| 221 | gbx%Ichan(1:gbx%nchan),gbx%surfem(1:gbx%nchan), & |
|---|
| 222 | gbx%co2,gbx%ch4,gbx%n2o,gbx%co, & |
|---|
| 223 | gbx_it) |
|---|
| 224 | ! --- Copy arrays without Npoints as dimension --- |
|---|
| 225 | gbx_it%dist_prmts_hydro = gbx%dist_prmts_hydro |
|---|
| 226 | gbx_it%dist_type_aero = gbx_it%dist_type_aero |
|---|
| 227 | call construct_cosp_vgrid(gbx_it,vgrid%Nlvgrid,vgrid%use_vgrid,vgrid%csat_vgrid,vgrid_it) |
|---|
| 228 | call construct_cosp_subgrid(Ni, Ncolumns, Nlevels, sgx_it) |
|---|
| 229 | call construct_cosp_sgradar(cfg,Ni,Ncolumns,Nlevels,N_HYDRO,sgradar_it) |
|---|
| 230 | call construct_cosp_sglidar(cfg,Ni,Ncolumns,Nlevels,N_HYDRO,PARASOL_NREFL,sglidar_it) |
|---|
| 231 | call construct_cosp_isccp(cfg,Ni,Ncolumns,Nlevels,isccp_it) |
|---|
| 232 | call construct_cosp_misr(cfg,Ni,misr_it) |
|---|
| 233 | call construct_cosp_radarstats(cfg,Ni,Ncolumns,vgrid%Nlvgrid,N_HYDRO,stradar_it) |
|---|
| 234 | call construct_cosp_lidarstats(cfg,Ni,Ncolumns,vgrid%Nlvgrid,N_HYDRO,PARASOL_NREFL,stlidar_it) |
|---|
| 235 | endif |
|---|
| 236 | ! --- Copy sections of arrays with Npoints as dimension --- |
|---|
| 237 | ix=(/i_first,i_last/) |
|---|
| 238 | iy=(/1,Ni/) |
|---|
| 239 | call cosp_gridbox_cpsection(ix,iy,gbx,gbx_it) |
|---|
| 240 | ! These serve as initialisation of *_it types |
|---|
| 241 | call cosp_subgrid_cpsection(ix,iy,sgx,sgx_it) |
|---|
| 242 | if (cfg%Lradar_sim) call cosp_sgradar_cpsection(ix,iy,sgradar,sgradar_it) |
|---|
| 243 | if (cfg%Llidar_sim) call cosp_sglidar_cpsection(ix,iy,sglidar,sglidar_it) |
|---|
| 244 | if (cfg%Lisccp_sim) call cosp_isccp_cpsection(ix,iy,isccp,isccp_it) |
|---|
| 245 | if (cfg%Lmisr_sim) call cosp_misr_cpsection(ix,iy,misr,misr_it) |
|---|
| 246 | if (cfg%Lradar_sim) call cosp_radarstats_cpsection(ix,iy,stradar,stradar_it) |
|---|
| 247 | if (cfg%Llidar_sim) call cosp_lidarstats_cpsection(ix,iy,stlidar,stlidar_it) |
|---|
| 248 | print *,'---------ix: ',ix |
|---|
| 249 | call cosp_iter(overlap,seed(ix(1):ix(2)),cfg,vgrid_it,gbx_it,sgx_it,sgradar_it, & |
|---|
| 250 | sglidar_it,isccp_it,misr_it,stradar_it,stlidar_it) |
|---|
| 251 | |
|---|
| 252 | ! --- Copy results to output structures --- |
|---|
| 253 | ! call cosp_gridbox_cphp(gbx_it,gbx) |
|---|
| 254 | ix=(/1,Ni/) |
|---|
| 255 | iy=(/i_first,i_last/) |
|---|
| 256 | call cosp_subgrid_cpsection(ix,iy,sgx_it,sgx) |
|---|
| 257 | if (cfg%Lradar_sim) call cosp_sgradar_cpsection(ix,iy,sgradar_it,sgradar) |
|---|
| 258 | if (cfg%Llidar_sim) call cosp_sglidar_cpsection(ix,iy,sglidar_it,sglidar) |
|---|
| 259 | if (cfg%Lisccp_sim) call cosp_isccp_cpsection(ix,iy,isccp_it,isccp) |
|---|
| 260 | if (cfg%Lmisr_sim) call cosp_misr_cpsection(ix,iy,misr_it,misr) |
|---|
| 261 | if (cfg%Lradar_sim) call cosp_radarstats_cpsection(ix,iy,stradar_it,stradar) |
|---|
| 262 | if (cfg%Llidar_sim) call cosp_lidarstats_cpsection(ix,iy,stlidar_it,stlidar) |
|---|
| 263 | enddo |
|---|
| 264 | ! Deallocate types |
|---|
| 265 | call free_cosp_gridbox(gbx_it,.true.) |
|---|
| 266 | call free_cosp_subgrid(sgx_it) |
|---|
| 267 | call free_cosp_vgrid(vgrid_it) |
|---|
| 268 | call free_cosp_sgradar(sgradar_it) |
|---|
| 269 | call free_cosp_sglidar(sglidar_it) |
|---|
| 270 | call free_cosp_isccp(isccp_it) |
|---|
| 271 | call free_cosp_misr(misr_it) |
|---|
| 272 | call free_cosp_radarstats(stradar_it) |
|---|
| 273 | call free_cosp_lidarstats(stlidar_it) |
|---|
| 274 | endif |
|---|
| 275 | deallocate(seed) |
|---|
| 276 | |
|---|
| 277 | |
|---|
| 278 | END SUBROUTINE COSP |
|---|
| 279 | |
|---|
| 280 | !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
|---|
| 281 | !--------------------- SUBROUTINE COSP_ITER ---------------------- |
|---|
| 282 | !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
|---|
| 283 | SUBROUTINE COSP_ITER(overlap,seed,cfg,vgrid,gbx,sgx,sgradar,sglidar,isccp,misr,stradar,stlidar) |
|---|
| 284 | |
|---|
| 285 | ! Arguments |
|---|
| 286 | integer,intent(in) :: overlap ! overlap type in SCOPS: 1=max, 2=rand, 3=max/rand |
|---|
| 287 | integer,dimension(:),intent(in) :: seed |
|---|
| 288 | type(cosp_config),intent(in) :: cfg ! Configuration options |
|---|
| 289 | type(cosp_vgrid),intent(in) :: vgrid ! Information on vertical grid of stats |
|---|
| 290 | type(cosp_gridbox),intent(inout) :: gbx |
|---|
| 291 | type(cosp_subgrid),intent(inout) :: sgx ! Subgrid info |
|---|
| 292 | type(cosp_sgradar),intent(inout) :: sgradar ! Output from radar simulator |
|---|
| 293 | type(cosp_sglidar),intent(inout) :: sglidar ! Output from lidar simulator |
|---|
| 294 | type(cosp_isccp),intent(inout) :: isccp ! Output from ISCCP simulator |
|---|
| 295 | type(cosp_misr),intent(inout) :: misr ! Output from MISR simulator |
|---|
| 296 | type(cosp_radarstats),intent(inout) :: stradar ! Summary statistics from radar simulator |
|---|
| 297 | type(cosp_lidarstats),intent(inout) :: stlidar ! Summary statistics from lidar simulator |
|---|
| 298 | |
|---|
| 299 | ! Local variables |
|---|
| 300 | integer :: Npoints ! Number of gridpoints |
|---|
| 301 | integer :: Ncolumns ! Number of subcolumns |
|---|
| 302 | integer :: Nlevels ! Number of levels |
|---|
| 303 | integer :: Nhydro ! Number of hydrometeors |
|---|
| 304 | integer :: Niter ! Number of calls to cosp_simulator |
|---|
| 305 | integer :: i,j,k |
|---|
| 306 | real,dimension(:,:),pointer :: column_frac_out ! Array with one column of frac_out |
|---|
| 307 | integer :: scops_debug=0 ! set to non-zero value to print out inputs for debugging in SCOPS |
|---|
| 308 | real,dimension(:, :),allocatable :: cca_scops,ls_p_rate,cv_p_rate, & |
|---|
| 309 | tca_scops ! Cloud cover in each model level (HORIZONTAL gridbox fraction) of total cloud. |
|---|
| 310 | ! Levels are from TOA to SURFACE. (nPoints, nLev) |
|---|
| 311 | real,dimension(:,:),allocatable :: frac_ls,prec_ls,frac_cv,prec_cv ! Cloud/Precipitation fraction in each model level |
|---|
| 312 | ! Levels are from SURFACE to TOA |
|---|
| 313 | type(cosp_sghydro) :: sghydro ! Subgrid info for hydrometeors en each iteration |
|---|
| 314 | |
|---|
| 315 | |
|---|
| 316 | !++++++++++ Dimensions ++++++++++++ |
|---|
| 317 | Npoints = gbx%Npoints |
|---|
| 318 | Ncolumns = gbx%Ncolumns |
|---|
| 319 | Nlevels = gbx%Nlevels |
|---|
| 320 | Nhydro = gbx%Nhydro |
|---|
| 321 | |
|---|
| 322 | |
|---|
| 323 | !++++++++++ Climate/NWP mode ++++++++++ |
|---|
| 324 | if (Ncolumns > 1) then |
|---|
| 325 | !++++++++++ Subgrid sampling ++++++++++ |
|---|
| 326 | ! Allocate arrays before calling SCOPS |
|---|
| 327 | allocate(frac_ls(Npoints,Nlevels),frac_cv(Npoints,Nlevels),prec_ls(Npoints,Nlevels),prec_cv(Npoints,Nlevels)) |
|---|
| 328 | allocate(tca_scops(Npoints,Nlevels),cca_scops(Npoints,Nlevels), & |
|---|
| 329 | ls_p_rate(Npoints,Nlevels),cv_p_rate(Npoints,Nlevels)) |
|---|
| 330 | ! Initialize to zero |
|---|
| 331 | frac_ls=0.0 |
|---|
| 332 | prec_ls=0.0 |
|---|
| 333 | frac_cv=0.0 |
|---|
| 334 | prec_cv=0.0 |
|---|
| 335 | ! Cloud fractions for SCOPS from TOA to SFC |
|---|
| 336 | tca_scops = gbx%tca(:,Nlevels:1:-1) |
|---|
| 337 | cca_scops = gbx%cca(:,Nlevels:1:-1) |
|---|
| 338 | |
|---|
| 339 | ! Call to SCOPS |
|---|
| 340 | ! strat and conv arrays are passed with levels from TOA to SURFACE. |
|---|
| 341 | call scops(Npoints,Nlevels,Ncolumns,seed,tca_scops,cca_scops,overlap,sgx%frac_out,scops_debug) |
|---|
| 342 | |
|---|
| 343 | ! temporarily use prec_ls/cv to transfer information about precipitation flux into prec_scops |
|---|
| 344 | if(gbx%use_precipitation_fluxes) then |
|---|
| 345 | ls_p_rate(:,Nlevels:1:-1)=gbx%rain_ls(:,1:Nlevels)+gbx%snow_ls(:,1:Nlevels)+gbx%grpl_ls(:,1:Nlevels) |
|---|
| 346 | cv_p_rate(:,Nlevels:1:-1)=gbx%rain_cv(:,1:Nlevels)+gbx%snow_cv(:,1:Nlevels) |
|---|
| 347 | else |
|---|
| 348 | ls_p_rate(:,Nlevels:1:-1)=gbx%mr_hydro(:,1:Nlevels,I_LSRAIN)+ & |
|---|
| 349 | gbx%mr_hydro(:,1:Nlevels,I_LSSNOW)+ & |
|---|
| 350 | gbx%mr_hydro(:,1:Nlevels,I_LSGRPL) |
|---|
| 351 | cv_p_rate(:,Nlevels:1:-1)=gbx%mr_hydro(:,1:Nlevels,I_CVRAIN)+ & |
|---|
| 352 | gbx%mr_hydro(:,1:Nlevels,I_CVSNOW) |
|---|
| 353 | endif |
|---|
| 354 | |
|---|
| 355 | call prec_scops(Npoints,Nlevels,Ncolumns,ls_p_rate,cv_p_rate,sgx%frac_out,sgx%prec_frac) |
|---|
| 356 | |
|---|
| 357 | ! Precipitation fraction |
|---|
| 358 | do j=1,Npoints,1 |
|---|
| 359 | do k=1,Nlevels,1 |
|---|
| 360 | do i=1,Ncolumns,1 |
|---|
| 361 | if (sgx%frac_out (j,i,Nlevels+1-k) .eq. 1) frac_ls(j,k)=frac_ls(j,k)+1. |
|---|
| 362 | if (sgx%frac_out (j,i,Nlevels+1-k) .eq. 2) frac_cv(j,k)=frac_cv(j,k)+1. |
|---|
| 363 | if (sgx%prec_frac(j,i,Nlevels+1-k) .eq. 1) prec_ls(j,k)=prec_ls(j,k)+1. |
|---|
| 364 | if (sgx%prec_frac(j,i,Nlevels+1-k) .eq. 2) prec_cv(j,k)=prec_cv(j,k)+1. |
|---|
| 365 | if (sgx%prec_frac(j,i,Nlevels+1-k) .eq. 3) then |
|---|
| 366 | prec_cv(j,k)=prec_cv(j,k)+1. |
|---|
| 367 | prec_ls(j,k)=prec_ls(j,k)+1. |
|---|
| 368 | endif |
|---|
| 369 | enddo !i |
|---|
| 370 | frac_ls(j,k)=frac_ls(j,k)/Ncolumns |
|---|
| 371 | frac_cv(j,k)=frac_cv(j,k)/Ncolumns |
|---|
| 372 | prec_ls(j,k)=prec_ls(j,k)/Ncolumns |
|---|
| 373 | prec_cv(j,k)=prec_cv(j,k)/Ncolumns |
|---|
| 374 | enddo !k |
|---|
| 375 | enddo !j |
|---|
| 376 | |
|---|
| 377 | ! Levels from SURFACE to TOA. |
|---|
| 378 | if (Npoints*Ncolumns*Nlevels < 10000) then |
|---|
| 379 | sgx%frac_out(1:Npoints,:,1:Nlevels) = sgx%frac_out(1:Npoints,:,Nlevels:1:-1) |
|---|
| 380 | sgx%prec_frac(1:Npoints,:,1:Nlevels) = sgx%prec_frac(1:Npoints,:,Nlevels:1:-1) |
|---|
| 381 | else |
|---|
| 382 | ! This is done within a loop (unvectorized) over nPoints to save memory |
|---|
| 383 | do j=1,Npoints |
|---|
| 384 | sgx%frac_out(j,:,1:Nlevels) = sgx%frac_out(j,:,Nlevels:1:-1) |
|---|
| 385 | sgx%prec_frac(j,:,1:Nlevels) = sgx%prec_frac(j,:,Nlevels:1:-1) |
|---|
| 386 | enddo |
|---|
| 387 | endif |
|---|
| 388 | |
|---|
| 389 | ! Deallocate arrays that will no longer be used |
|---|
| 390 | deallocate(tca_scops,cca_scops,ls_p_rate,cv_p_rate) |
|---|
| 391 | |
|---|
| 392 | ! Populate the subgrid arrays |
|---|
| 393 | call construct_cosp_sghydro(Npoints,Ncolumns,Nlevels,Nhydro,sghydro) |
|---|
| 394 | do k=1,Ncolumns |
|---|
| 395 | !--------- Mixing ratios for clouds and Reff for Clouds and precip ------- |
|---|
| 396 | column_frac_out => sgx%frac_out(:,k,:) |
|---|
| 397 | where (column_frac_out == 1) !+++++++++++ LS clouds ++++++++ |
|---|
| 398 | sghydro%mr_hydro(:,k,:,I_LSCLIQ) = gbx%mr_hydro(:,:,I_LSCLIQ) |
|---|
| 399 | sghydro%mr_hydro(:,k,:,I_LSCICE) = gbx%mr_hydro(:,:,I_LSCICE) |
|---|
| 400 | |
|---|
| 401 | sghydro%Reff(:,k,:,I_LSCLIQ) = gbx%Reff(:,:,I_LSCLIQ) |
|---|
| 402 | sghydro%Reff(:,k,:,I_LSCICE) = gbx%Reff(:,:,I_LSCICE) |
|---|
| 403 | sghydro%Reff(:,k,:,I_LSRAIN) = gbx%Reff(:,:,I_LSRAIN) |
|---|
| 404 | sghydro%Reff(:,k,:,I_LSSNOW) = gbx%Reff(:,:,I_LSSNOW) |
|---|
| 405 | sghydro%Reff(:,k,:,I_LSGRPL) = gbx%Reff(:,:,I_LSGRPL) |
|---|
| 406 | elsewhere (column_frac_out == 2) !+++++++++++ CONV clouds ++++++++ |
|---|
| 407 | sghydro%mr_hydro(:,k,:,I_CVCLIQ) = gbx%mr_hydro(:,:,I_CVCLIQ) |
|---|
| 408 | sghydro%mr_hydro(:,k,:,I_CVCICE) = gbx%mr_hydro(:,:,I_CVCICE) |
|---|
| 409 | |
|---|
| 410 | sghydro%Reff(:,k,:,I_CVCLIQ) = gbx%Reff(:,:,I_CVCLIQ) |
|---|
| 411 | sghydro%Reff(:,k,:,I_CVCICE) = gbx%Reff(:,:,I_CVCICE) |
|---|
| 412 | sghydro%Reff(:,k,:,I_CVRAIN) = gbx%Reff(:,:,I_CVRAIN) |
|---|
| 413 | sghydro%Reff(:,k,:,I_CVSNOW) = gbx%Reff(:,:,I_CVSNOW) |
|---|
| 414 | end where |
|---|
| 415 | !--------- Precip ------- |
|---|
| 416 | if (.not. gbx%use_precipitation_fluxes) then |
|---|
| 417 | where (column_frac_out == 1) !+++++++++++ LS Precipitation ++++++++ |
|---|
| 418 | sghydro%mr_hydro(:,k,:,I_LSRAIN) = gbx%mr_hydro(:,:,I_LSRAIN) |
|---|
| 419 | sghydro%mr_hydro(:,k,:,I_LSSNOW) = gbx%mr_hydro(:,:,I_LSSNOW) |
|---|
| 420 | sghydro%mr_hydro(:,k,:,I_LSGRPL) = gbx%mr_hydro(:,:,I_LSGRPL) |
|---|
| 421 | elsewhere (column_frac_out == 2) !+++++++++++ CONV Precipitation ++++++++ |
|---|
| 422 | sghydro%mr_hydro(:,k,:,I_CVRAIN) = gbx%mr_hydro(:,:,I_CVRAIN) |
|---|
| 423 | sghydro%mr_hydro(:,k,:,I_CVSNOW) = gbx%mr_hydro(:,:,I_CVSNOW) |
|---|
| 424 | end where |
|---|
| 425 | endif |
|---|
| 426 | enddo |
|---|
| 427 | ! convert the mixing ratio and precipitation flux from gridbox mean to the fraction-based values |
|---|
| 428 | do k=1,Nlevels |
|---|
| 429 | do j=1,Npoints |
|---|
| 430 | !--------- Clouds ------- |
|---|
| 431 | if (frac_ls(j,k) .ne. 0.) then |
|---|
| 432 | sghydro%mr_hydro(j,:,k,I_LSCLIQ) = sghydro%mr_hydro(j,:,k,I_LSCLIQ)/frac_ls(j,k) |
|---|
| 433 | sghydro%mr_hydro(j,:,k,I_LSCICE) = sghydro%mr_hydro(j,:,k,I_LSCICE)/frac_ls(j,k) |
|---|
| 434 | endif |
|---|
| 435 | if (frac_cv(j,k) .ne. 0.) then |
|---|
| 436 | sghydro%mr_hydro(j,:,k,I_CVCLIQ) = sghydro%mr_hydro(j,:,k,I_CVCLIQ)/frac_cv(j,k) |
|---|
| 437 | sghydro%mr_hydro(j,:,k,I_CVCICE) = sghydro%mr_hydro(j,:,k,I_CVCICE)/frac_cv(j,k) |
|---|
| 438 | endif |
|---|
| 439 | !--------- Precip ------- |
|---|
| 440 | if (gbx%use_precipitation_fluxes) then |
|---|
| 441 | if (prec_ls(j,k) .ne. 0.) then |
|---|
| 442 | gbx%rain_ls(j,k) = gbx%rain_ls(j,k)/prec_ls(j,k) |
|---|
| 443 | gbx%snow_ls(j,k) = gbx%snow_ls(j,k)/prec_ls(j,k) |
|---|
| 444 | gbx%grpl_ls(j,k) = gbx%grpl_ls(j,k)/prec_ls(j,k) |
|---|
| 445 | endif |
|---|
| 446 | if (prec_cv(j,k) .ne. 0.) then |
|---|
| 447 | gbx%rain_cv(j,k) = gbx%rain_cv(j,k)/prec_cv(j,k) |
|---|
| 448 | gbx%snow_cv(j,k) = gbx%snow_cv(j,k)/prec_cv(j,k) |
|---|
| 449 | endif |
|---|
| 450 | else |
|---|
| 451 | if (prec_ls(j,k) .ne. 0.) then |
|---|
| 452 | sghydro%mr_hydro(j,:,k,I_LSRAIN) = sghydro%mr_hydro(j,:,k,I_LSRAIN)/prec_ls(j,k) |
|---|
| 453 | sghydro%mr_hydro(j,:,k,I_LSSNOW) = sghydro%mr_hydro(j,:,k,I_LSSNOW)/prec_ls(j,k) |
|---|
| 454 | sghydro%mr_hydro(j,:,k,I_LSGRPL) = sghydro%mr_hydro(j,:,k,I_LSGRPL)/prec_ls(j,k) |
|---|
| 455 | endif |
|---|
| 456 | if (prec_cv(j,k) .ne. 0.) then |
|---|
| 457 | sghydro%mr_hydro(j,:,k,I_CVRAIN) = sghydro%mr_hydro(j,:,k,I_CVRAIN)/prec_cv(j,k) |
|---|
| 458 | sghydro%mr_hydro(j,:,k,I_CVSNOW) = sghydro%mr_hydro(j,:,k,I_CVSNOW)/prec_cv(j,k) |
|---|
| 459 | endif |
|---|
| 460 | endif |
|---|
| 461 | enddo !k |
|---|
| 462 | enddo !j |
|---|
| 463 | deallocate(frac_ls,prec_ls,frac_cv,prec_cv) |
|---|
| 464 | |
|---|
| 465 | if (gbx%use_precipitation_fluxes) then |
|---|
| 466 | ! convert precipitation flux into mixing ratio |
|---|
| 467 | call pf_to_mr(Npoints,Nlevels,Ncolumns,gbx%rain_ls,gbx%snow_ls,gbx%grpl_ls, & |
|---|
| 468 | gbx%rain_cv,gbx%snow_cv,sgx%prec_frac,gbx%p,gbx%T, & |
|---|
| 469 | sghydro%mr_hydro(:,:,:,I_LSRAIN),sghydro%mr_hydro(:,:,:,I_LSSNOW),sghydro%mr_hydro(:,:,:,I_LSGRPL), & |
|---|
| 470 | sghydro%mr_hydro(:,:,:,I_CVRAIN),sghydro%mr_hydro(:,:,:,I_CVSNOW)) |
|---|
| 471 | endif |
|---|
| 472 | !++++++++++ CRM mode ++++++++++ |
|---|
| 473 | else |
|---|
| 474 | sghydro%mr_hydro(:,1,:,:) = gbx%mr_hydro |
|---|
| 475 | sghydro%Reff(:,1,:,:) = gbx%Reff |
|---|
| 476 | !--------- Clouds ------- |
|---|
| 477 | where ((gbx%dtau_s > 0.0)) |
|---|
| 478 | sgx%frac_out(:,1,:) = 1 ! Subgrid cloud array. Dimensions (Npoints,Ncolumns,Nlevels) |
|---|
| 479 | endwhere |
|---|
| 480 | endif ! Ncolumns > 1 |
|---|
| 481 | |
|---|
| 482 | |
|---|
| 483 | !++++++++++ Simulator ++++++++++ |
|---|
| 484 | call cosp_simulator(gbx,sgx,sghydro,cfg,vgrid,sgradar,sglidar,isccp,misr,stradar,stlidar) |
|---|
| 485 | |
|---|
| 486 | ! Deallocate subgrid arrays |
|---|
| 487 | call free_cosp_sghydro(sghydro) |
|---|
| 488 | END SUBROUTINE COSP_ITER |
|---|
| 489 | |
|---|
| 490 | END MODULE MOD_COSP |
|---|
