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