source: LMDZ5/branches/IPSLCM5A2.1_ISO/libf/phyiso/cosp/phys_cosp.F90 @ 5360

Last change on this file since 5360 was 3331, checked in by acozic, 7 years ago

Add modification for isotopes

  • Property svn:executable set to *
File size: 18.7 KB
Line 
1! Simulateur COSP : Cfmip Observation Simulator Package
2
3! ISCCP, Radar (QuickBeam), Lidar et Parasol (ACTSIM), MISR, RTTOVS
4!Idelkadi Abderrahmane Aout-Septembre 2009 First Version
5!Idelkadi Abderrahmane Nov 2015 version v1.4.0
6
7  subroutine phys_cosp( itap,dtime,freq_cosp, &
8                        ok_mensuelCOSP,ok_journeCOSP,ok_hfCOSP, &
9                        ecrit_mth,ecrit_day,ecrit_hf, ok_all_xml, &
10                        Nptslmdz,Nlevlmdz,lon,lat, presnivs,overlaplmdz,sunlit, &
11                        ref_liq,ref_ice,fracTerLic,u_wind,v_wind,phis,phi,ph,p,skt,t, &
12                        sh,rh,tca,cca,mr_lsliq,mr_lsice,fl_lsrainI,fl_lssnowI, &
13                        fl_ccrainI,fl_ccsnowI,mr_ozone,dtau_s,dem_s)
14
15!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
16!!!! Inputs :
17! itap,                                 !Increment de la physiq
18! dtime,                                !Pas de temps physiq
19! overlap,                              !Overlap type in SCOPS
20! Npoints,                              !Nb de points de la grille physiq
21! Nlevels,                              !Nb de niveaux verticaux
22! Ncolumns,                             !Number of subcolumns
23! lon,lat,                              !Longitudes et latitudes de la grille LMDZ
24! ref_liq,ref_ice,                      !Rayons effectifs des particules liq et ice (en microm)
25! fracTerLic,                           !Fraction terre a convertir en masque
26! u_wind,v_wind,                        !Vents a 10m ???
27! phi,                                  !Geopotentiel
28! phis,                                 !Geopotentiel sol
29! ph,                                   !pression pour chaque inter-couche
30! p,                                    !Pression aux milieux des couches
31! skt,t,                                !Temp au sol et temp 3D
32! sh,                                   !Humidite specifique
33! rh,                                   !Humidite relatif
34! tca,                                  !Fraction nuageuse
35! cca                                   !Fraction nuageuse convective
36! mr_lsliq,                             !Liq Cloud water content
37! mr_lsice,                             !Ice Cloud water content
38! mr_ccliq,                             !Convective Cloud Liquid water content 
39! mr_ccice,                             !Cloud ice water content
40! fl_lsrain,                            !Large scale precipitation lic
41! fl_lssnow,                            !Large scale precipitation ice
42! fl_ccrain,                            !Convective precipitation lic
43! fl_ccsnow,                            !Convective precipitation ice
44! mr_ozone,                             !Concentration ozone (Kg/Kg)
45! dem_s                                 !Cloud optical emissivity
46! dtau_s                                !Cloud optical thickness
47! emsfc_lw = 1.                         !Surface emissivity dans radlwsw.F90
48
49!!! Outputs :
50! calipso2D,                            !Lidar Low/heigh/Mean/Total-level Cloud Fraction
51! calipso3D,                            !Lidar Cloud Fraction (532 nm)
52! cfadlidar,                            !Lidar Scattering Ratio CFAD (532 nm)
53! parasolrefl,                          !PARASOL-like mono-directional reflectance
54! atb,                                  !Lidar Attenuated Total Backscatter (532 nm)
55! betamol,                              !Lidar Molecular Backscatter (532 nm)
56! cfaddbze,                             !Radar Reflectivity Factor CFAD (94 GHz)
57! clcalipso2,                           !Cloud frequency of occurrence as seen by CALIPSO but not CloudSat
58! dbze,                                 !Efective_reflectivity_factor
59! cltlidarradar,                        !Lidar and Radar Total Cloud Fraction
60! clMISR,                               !Cloud Fraction as Calculated by the MISR Simulator
61! clisccp2,                             !Cloud Fraction as Calculated by the ISCCP Simulator
62! boxtauisccp,                          !Optical Depth in Each Column as Calculated by the ISCCP Simulator
63! boxptopisccp,                         !Cloud Top Pressure in Each Column as Calculated by the ISCCP Simulator
64! tclisccp,                             !Total Cloud Fraction as Calculated by the ISCCP Simulator
65! ctpisccp,                             !Mean Cloud Top Pressure as Calculated by the ISCCP Simulator
66! tauisccp,                             !Mean Optical Depth as Calculated by the ISCCP Simulator
67! albisccp,                             !Mean Cloud Albedo as Calculated by the ISCCP Simulator
68! meantbisccp,                          !Mean all-sky 10.5 micron brightness temperature as calculated by the ISCCP Simulator
69! meantbclrisccp                        !Mean clear-sky 10.5 micron brightness temperature as calculated by the ISCCP Simulator
70
71!!! AI rajouter les nouvelles sorties
72!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
73
74!! AI rajouter
75#include "cosp_defs.h"
76  USE MOD_COSP_CONSTANTS
77  USE MOD_COSP_TYPES
78  USE MOD_COSP
79  USE mod_phys_lmdz_para
80  USE mod_grid_phy_lmdz
81  use ioipsl
82  use iophy
83  use cosp_output_mod
84  use cosp_output_write_mod
85!  use MOD_COSP_Modis_Simulator, only : cosp_modis
86
87  IMPLICIT NONE
88
89  ! Local variables
90  character(len=64),PARAMETER  :: cosp_input_nl='cosp_input_nl.txt'
91  character(len=64),PARAMETER  :: cosp_output_nl='cosp_output_nl.txt'
92  integer, save :: isccp_topheight,isccp_topheight_direction,overlap
93  integer,save  :: Ncolumns     ! Number of subcolumns in SCOPS
94  integer, save :: Npoints      ! Number of gridpoints
95!$OMP THREADPRIVATE(Npoints)
96  integer, save :: Nlevels      ! Number of levels
97  Integer :: Nptslmdz,Nlevlmdz ! Nb de points issus de physiq.F
98  integer, save :: Nlr          ! Number of levels in statistical outputs
99  integer, save :: Npoints_it   ! Max number of gridpoints to be processed in one iteration
100  integer :: i
101  type(cosp_config),save :: cfg   ! Configuration options
102!$OMP THREADPRIVATE(cfg)
103  type(cosp_gridbox) :: gbx ! Gridbox information. Input for COSP
104  type(cosp_subgrid) :: sgx     ! Subgrid outputs
105  type(cosp_sgradar) :: sgradar ! Output from radar simulator
106  type(cosp_sglidar) :: sglidar ! Output from lidar simulator
107  type(cosp_isccp)   :: isccp   ! Output from ISCCP simulator
108!! AI rajout modis
109  type(cosp_modis)   :: modis   ! Output from MODIS simulator
110!!
111  type(cosp_misr)    :: misr    ! Output from MISR simulator
112!! AI rajout rttovs
113!  type(cosp_rttov)   :: rttov   ! Output from RTTOV
114!!
115  type(cosp_vgrid)   :: vgrid   ! Information on vertical grid of stats
116  type(cosp_radarstats) :: stradar ! Summary statistics from radar simulator
117  type(cosp_lidarstats) :: stlidar ! Summary statistics from lidar simulator
118
119  integer :: t0,t1,count_rate,count_max
120  integer :: Nlon,Nlat
121  real,save :: radar_freq,k2,ZenAng,co2,ch4,n2o,co,emsfc_lw
122!$OMP THREADPRIVATE(emsfc_lw)
123  integer,dimension(RTTOV_MAX_CHANNELS),save :: Channels
124  real,dimension(RTTOV_MAX_CHANNELS),save :: Surfem
125  integer, save :: surface_radar,use_mie_tables,use_gas_abs,do_ray,melt_lay
126  integer, save :: Nprmts_max_hydro,Naero,Nprmts_max_aero,lidar_ice_type
127  integer, save :: platform,satellite,Instrument,Nchannels
128  logical, save :: use_vgrid,csat_vgrid,use_precipitation_fluxes,use_reff
129
130! Declaration necessaires pour les sorties IOIPSL
131  integer :: ii
132  real    :: ecrit_day,ecrit_hf,ecrit_mth
133  logical :: ok_mensuelCOSP,ok_journeCOSP,ok_hfCOSP, ok_all_xml
134
135  logical, save :: debut_cosp=.true.
136!$OMP THREADPRIVATE(debut_cosp)
137
138!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! Input variables from LMDZ-GCM
139  integer                         :: overlaplmdz   !  overlap type: 1=max, 2=rand, 3=max/rand ! cosp input (output lmdz)
140  real,dimension(Nptslmdz,Nlevlmdz) :: height,phi,p,ph,T,sh,rh,tca,cca,mr_lsliq,mr_lsice,mr_ccliq,mr_ccice, &
141                                     fl_lsrain,fl_lssnow,fl_ccrain,fl_ccsnow,fl_lsgrpl, &
142                                     zlev,zlev_half,mr_ozone,radliq,radice,dtau_s,dem_s,ref_liq,ref_ice
143  real,dimension(Nptslmdz,Nlevlmdz) ::  fl_lsrainI,fl_lssnowI,fl_ccrainI,fl_ccsnowI
144  real,dimension(Nptslmdz)        :: lon,lat,skt,fracTerLic,u_wind,v_wind,phis,sunlit         
145  real,dimension(Nlevlmdz)        :: presnivs
146  integer                         :: itap,k,ip
147  real                            :: dtime,freq_cosp
148  real,dimension(2)               :: time_bnds
149 
150   namelist/COSP_INPUT/overlap,isccp_topheight,isccp_topheight_direction, &
151              npoints_it,ncolumns,use_vgrid,nlr,csat_vgrid, &
152              radar_freq,surface_radar,use_mie_tables, &
153              use_gas_abs,do_ray,melt_lay,k2,Nprmts_max_hydro,Naero,Nprmts_max_aero, &
154              lidar_ice_type,use_precipitation_fluxes,use_reff, &
155              platform,satellite,Instrument,Nchannels, &
156              Channels,Surfem,ZenAng,co2,ch4,n2o,co
157
158!---------------- End of declaration of variables --------------
159   
160
161!+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
162! Read namelist with COSP inputs
163!+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
164
165 if (debut_cosp) then
166  NPoints=Nptslmdz
167  Nlevels=Nlevlmdz
168 
169! Lecture du namelist input
170  CALL read_cosp_input
171
172! Clefs Outputs
173  call read_cosp_output_nl(cosp_output_nl,cfg)
174
175    if (overlaplmdz.ne.overlap) then
176       print*,'Attention overlaplmdz different de overlap lu dans namelist '
177    endif
178   print*,'Fin lecture Namelists, debut_cosp =',debut_cosp
179
180  print*,' Cles des differents simulateurs cosp :'
181  print*,'Lradar_sim,Llidar_sim,Lisccp_sim,Lmisr_sim,Lmodis_sim,Lrttov_sim', &
182          cfg%Lradar_sim,cfg%Llidar_sim,cfg%Lisccp_sim,cfg%Lmisr_sim,cfg%Lmodis_sim,cfg%Lrttov_sim
183
184  endif ! debut_cosp
185
186  time_bnds(1) = dtime-dtime/2.
187  time_bnds(2) = dtime+dtime/2.
188
189!  print*,'Debut phys_cosp itap,dtime,freq_cosp,ecrit_mth,ecrit_day,ecrit_hf ', &
190!          itap,dtime,freq_cosp,ecrit_mth,ecrit_day,ecrit_hf
191!+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
192! Allocate memory for gridbox type
193!+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
194        print *, 'Allocating memory for gridbox type...'
195!! AI
196!        call construct_cosp_gridbox(dble(itap),radar_freq,surface_radar,use_mie_tables,use_gas_abs,do_ray,melt_lay,k2, &
197!                                    Npoints,Nlevels,Ncolumns,N_HYDRO,Nprmts_max_hydro,Naero,Nprmts_max_aero,Npoints_it, &
198!                                    lidar_ice_type,isccp_topheight,isccp_topheight_direction,overlap,emsfc_lw, &
199!                                    use_precipitation_fluxes,use_reff, &
200!                                    Platform,Satellite,Instrument,Nchannels,ZenAng, &
201!                                    channels(1:Nchannels),surfem(1:Nchannels),co2,ch4,n2o,co,gbx)
202! Surafce emissivity
203        emsfc_lw = 1.
204
205        call construct_cosp_gridbox(dtime,time_bnds,radar_freq,surface_radar,use_mie_tables,use_gas_abs, &
206                                    do_ray,melt_lay,k2, &
207                                    Npoints,Nlevels,Ncolumns,N_HYDRO,Nprmts_max_hydro,Naero,Nprmts_max_aero,Npoints_it, &
208                                    lidar_ice_type,isccp_topheight,isccp_topheight_direction,overlap,emsfc_lw, &
209                                    use_precipitation_fluxes,use_reff, &
210                                    Platform,Satellite,Instrument,Nchannels,ZenAng, &
211                                    channels(1:Nchannels),surfem(1:Nchannels),co2,ch4,n2o,co,gbx)
212       
213!+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
214! Here code to populate input structure
215!+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
216
217        print *, 'Populating input structure...'
218        gbx%longitude = lon
219        gbx%latitude = lat
220
221        gbx%p = p !
222        gbx%ph = ph
223        gbx%zlev = phi/9.81
224
225        zlev_half(:,1) = phis(:)/9.81
226        do k = 2, Nlevels
227          do ip = 1, Npoints
228           zlev_half(ip,k) = phi(ip,k)/9.81 + &
229               (phi(ip,k)-phi(ip,k-1))/9.81 * (ph(ip,k)-p(ip,k)) / (p(ip,k)-p(ip,k-1))
230          enddo
231        enddo
232        gbx%zlev_half = zlev_half
233
234        gbx%T = T
235        gbx%q = rh*100.
236        gbx%sh = sh
237! On ne veut pas que cosp distingue les nuages stratiformes et convectifs
238! on passe les contenus totaux (conv+strat)
239        gbx%cca = 0. !convective_cloud_amount (1)
240        gbx%tca = tca ! total_cloud_amount (1)
241        gbx%psfc = ph(:,1) !pression de surface
242        gbx%skt  = skt !Skin temperature (K)
243
244        do ip = 1, Npoints
245          if (fracTerLic(ip).ge.0.5) then
246             gbx%land(ip) = 1.
247          else
248             gbx%land(ip) = 0.
249          endif
250        enddo
251        gbx%mr_ozone  = mr_ozone !mass_fraction_of_ozone_in_air (kg/kg)
252! A voir l equivalent LMDZ (u10m et v10m)
253        gbx%u_wind  = u_wind !eastward_wind (m s-1)
254        gbx%v_wind  = v_wind !northward_wind
255
256! sunlit calcule a partir de la fraction d ensoleillement par jour
257!      do ip = 1, Npoints
258!        if (sunlit(ip).le.0.) then
259!           gbx%sunlit(ip)=0.
260!        else
261!           gbx%sunlit(ip)=1.
262!        endif
263!      enddo
264       gbx%sunlit=sunlit
265
266! A voir l equivalent LMDZ
267  mr_ccliq = 0.0
268  mr_ccice = 0.0
269        gbx%mr_hydro(:,:,I_LSCLIQ) = mr_lsliq !mixing_ratio_large_scale_cloud_liquid (kg/kg)
270        gbx%mr_hydro(:,:,I_LSCICE) = mr_lsice !mixing_ratio_large_scale_cloud_ic
271        gbx%mr_hydro(:,:,I_CVCLIQ) = mr_ccliq !mixing_ratio_convective_cloud_liquid
272        gbx%mr_hydro(:,:,I_CVCICE) = mr_ccice !mixing_ratio_convective_cloud_ice
273! A revoir
274        fl_lsrain = fl_lsrainI + fl_ccrainI
275        fl_lssnow = fl_lssnowI + fl_ccsnowI
276        gbx%rain_ls = fl_lsrain !flux_large_scale_cloud_rain (kg m^-2 s^-1)
277        gbx%snow_ls = fl_lssnow !flux_large_scale_cloud_snow
278!  A voir l equivalent LMDZ
279        fl_lsgrpl=0.
280        fl_ccsnow = 0.
281        fl_ccrain = 0.
282        gbx%grpl_ls = fl_lsgrpl  !flux_large_scale_cloud_graupel
283        gbx%rain_cv = fl_ccrain  !flux_convective_cloud_rain
284        gbx%snow_cv = fl_ccsnow  !flux_convective_cloud_snow
285
286     gbx%Reff(:,:,I_LSCLIQ) = ref_liq*1e-6
287     gbx%Reff(:,:,I_LSCICE) = ref_ice*1e-6
288!! AI A revoir
289     gbx%Reff(:,:,I_CVCLIQ) = ref_liq*1e-6
290     gbx%Reff(:,:,I_CVCICE) = ref_ice*1e-6
291
292        ! ISCCP simulator
293        gbx%dtau_s   = dtau_s
294        gbx%dtau_c   = 0.
295        gbx%dem_s    = dem_s
296        gbx%dem_c    = 0.
297
298!+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
299        ! Define new vertical grid
300!+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
301        print *, 'Defining new vertical grid...'
302        call construct_cosp_vgrid(gbx,Nlr,use_vgrid,csat_vgrid,vgrid)
303
304!+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
305       ! Allocate memory for other types
306!+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
307        print *, 'Allocating memory for other types...'
308        call construct_cosp_subgrid(Npoints, Ncolumns, Nlevels, sgx)
309        call construct_cosp_sgradar(cfg,Npoints,Ncolumns,Nlevels,N_HYDRO,sgradar)
310        call construct_cosp_radarstats(cfg,Npoints,Ncolumns,vgrid%Nlvgrid,N_HYDRO,stradar)
311        call construct_cosp_sglidar(cfg,Npoints,Ncolumns,Nlevels,N_HYDRO,PARASOL_NREFL,sglidar)
312        call construct_cosp_lidarstats(cfg,Npoints,Ncolumns,vgrid%Nlvgrid,N_HYDRO,PARASOL_NREFL,stlidar)
313        call construct_cosp_isccp(cfg,Npoints,Ncolumns,Nlevels,isccp)
314!! AI rajout
315        call construct_cosp_modis(cfg,Npoints,modis)
316!!
317        call construct_cosp_misr(cfg,Npoints,misr)
318!        call construct_cosp_rttov(cfg,Npoints,Nchannels,rttov)
319
320!+++++++++++++ Open output files and define output files axis !+++++++++++++
321     if (debut_cosp) then
322
323      !$OMP MASTER
324        print *, ' Open outpts files and define axis'
325        call cosp_output_open(Nlevlmdz, Ncolumns, presnivs, dtime, freq_cosp, &
326                              ok_mensuelCOSP, ok_journeCOSP, ok_hfCOSP, ok_all_xml, &
327                              ecrit_mth, ecrit_day, ecrit_hf, use_vgrid, vgrid)
328      !$OMP END MASTER
329      !$OMP BARRIER
330        debut_cosp=.false.
331      endif ! debut_cosp
332!+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
333        ! Call simulator
334!+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
335        print *, 'Calling simulator...'
336!! AI
337!        call cosp(overlap,Ncolumns,cfg,vgrid,gbx,sgx,sgradar,sglidar,isccp,misr,stradar,stlidar)
338!#ifdef RTTOV
339!        call cosp(overlap,Ncolumns,cfg,vgrid,gbx,sgx,sgradar,sglidar,isccp,misr,modis,rttov,stradar,stlidar)
340!#else
341        call cosp(overlap,Ncolumns,cfg,vgrid,gbx,sgx,sgradar,sglidar,isccp,misr,modis,stradar,stlidar)
342!#endif
343!!
344
345!+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
346
347!!!!!!!!!!!!!!!!!! Ecreture des sorties Cosp !!!!!!!!!!!!!!r!!!!!!:!!!!!
348
349       print *, 'Calling write output'
350        call cosp_output_write(Nlevlmdz, Npoints, Ncolumns, itap, dtime, freq_COSP, &
351                               cfg, gbx, vgrid, sglidar, sgradar, stlidar, stradar, &
352                               isccp, misr, modis)
353
354!+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
355        ! Deallocate memory in derived types
356!+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
357        print *, 'Deallocating memory...'
358        call free_cosp_gridbox(gbx)
359        call free_cosp_subgrid(sgx)
360        call free_cosp_sgradar(sgradar)
361        call free_cosp_radarstats(stradar)
362        call free_cosp_sglidar(sglidar)
363        call free_cosp_lidarstats(stlidar)
364        call free_cosp_isccp(isccp)
365        call free_cosp_misr(misr)
366!! AI
367        call free_cosp_modis(modis)
368!        call free_cosp_rttov(rttov)
369!!
370        call free_cosp_vgrid(vgrid) 
371 
372!+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
373  ! Time in s. Only for testing purposes
374!  call system_clock(t1,count_rate,count_max)
375!  print *,(t1-t0)*1.0/count_rate
376 
377  CONTAINS
378 
379  SUBROUTINE read_cosp_input
380   
381    IF (is_master) THEN
382      OPEN(10,file=cosp_input_nl,status='old')
383      READ(10,nml=cosp_input)
384      CLOSE(10)
385    ENDIF
386    CALL bcast(overlap)
387    CALL bcast(isccp_topheight)
388    CALL bcast(isccp_topheight_direction)
389    CALL bcast(npoints_it)
390    CALL bcast(ncolumns)
391    CALL bcast(use_vgrid)
392    CALL bcast(nlr)
393    CALL bcast(csat_vgrid)
394    CALL bcast(radar_freq)
395    CALL bcast(surface_radar)
396    CALL bcast(use_mie_tables)
397    CALL bcast(use_gas_abs)
398    CALL bcast(do_ray)
399    CALL bcast(melt_lay)
400    CALL bcast(k2)
401    CALL bcast(Nprmts_max_hydro)
402    CALL bcast(Naero)
403    CALL bcast(Nprmts_max_aero)
404    CALL bcast(lidar_ice_type)
405    CALL bcast(use_precipitation_fluxes)
406    CALL bcast(use_reff)
407    CALL bcast(platform)
408    CALL bcast(satellite)
409    CALL bcast(Instrument)
410    CALL bcast(Nchannels)
411    CALL bcast(Channels)
412    CALL bcast(Surfem)
413    CALL bcast(ZenAng)
414    CALL bcast(co2)
415    CALL bcast(ch4)
416    CALL bcast(n2o)
417    CALL bcast(co)
418!$OMP BARRIER 
419  END SUBROUTINE read_cosp_input
420
421end subroutine phys_cosp
Note: See TracBrowser for help on using the repository browser.