1 | ! |
---|
2 | ! $Id$ |
---|
3 | ! |
---|
4 | MODULE etat0phys |
---|
5 | ! |
---|
6 | !******************************************************************************* |
---|
7 | ! Purpose: Create physical initial state using atmospheric fields from a |
---|
8 | ! database of atmospheric to initialize the model. |
---|
9 | !------------------------------------------------------------------------------- |
---|
10 | ! Comments: |
---|
11 | ! |
---|
12 | ! * This module is designed to work for Earth (and with ioipsl) |
---|
13 | ! |
---|
14 | ! * etat0phys_netcdf routine can access to NetCDF data through subroutines: |
---|
15 | ! "start_init_phys" for variables contained in file "ECPHY.nc": |
---|
16 | ! 'ST' : Surface temperature |
---|
17 | ! 'CDSW' : Soil moisture |
---|
18 | ! "start_init_orog" for variables contained in file "Relief.nc": |
---|
19 | ! 'RELIEF' : High resolution orography |
---|
20 | ! |
---|
21 | ! * The land mask and corresponding weights can be: |
---|
22 | ! 1) computed using the ocean mask from the ocean model (to ensure ocean |
---|
23 | ! fractions are the same for atmosphere and ocean) for coupled runs. |
---|
24 | ! File name: "o2a.nc" ; variable name: "OceMask" |
---|
25 | ! 2) computed from topography file "Relief.nc" for forced runs. |
---|
26 | ! |
---|
27 | ! * Allowed values for read_climoz flag are 0, 1 and 2: |
---|
28 | ! 0: do not read an ozone climatology |
---|
29 | ! 1: read a single ozone climatology that will be used day and night |
---|
30 | ! 2: read two ozone climatologies, the average day and night climatology |
---|
31 | ! and the daylight climatology |
---|
32 | !------------------------------------------------------------------------------- |
---|
33 | ! * There is a big mess with the longitude size. Should it be iml or iml+1 ? |
---|
34 | ! I have chosen to use the iml+1 as an argument to this routine and we declare |
---|
35 | ! internaly smaller fields when needed. This needs to be cleared once and for |
---|
36 | ! all in LMDZ. A convention is required. |
---|
37 | !------------------------------------------------------------------------------- |
---|
38 | |
---|
39 | USE ioipsl, ONLY: flininfo, flinopen, flinget, flinclo |
---|
40 | USE assert_eq_m, ONLY: assert_eq |
---|
41 | USE dimphy, ONLY: klon |
---|
42 | USE conf_dat_m, ONLY: conf_dat2d |
---|
43 | USE phys_state_var_mod, ONLY: zmea, zstd, zsig, zgam, zthe, zpic, zval, z0m, & |
---|
44 | solsw, solswfdiff, radsol, t_ancien, wake_deltat, wake_s, rain_fall, qsol, z0h, & |
---|
45 | sollw,sollwdown, rugoro, q_ancien, wake_deltaq, wake_pe, snow_fall, ratqs,w01, & |
---|
46 | sig1, ftsol, clwcon, fm_therm, wake_Cstar, pctsrf, entr_therm,radpas, f0,& |
---|
47 | zmax0,fevap, rnebcon,falb_dir, falb_dif, wake_fip, agesno, detr_therm, pbl_tke, & |
---|
48 | phys_state_var_init, ql_ancien, qs_ancien, prlw_ancien, prsw_ancien, & |
---|
49 | prw_ancien, u10m,v10m, treedrg, u_ancien, v_ancien, wake_delta_pbl_TKE, wake_dens, & |
---|
50 | ale_bl, ale_bl_trig, alp_bl, & |
---|
51 | ale_wake, ale_bl_stat, AWAKE_S, & |
---|
52 | cf_ancien, rvc_ancien |
---|
53 | |
---|
54 | USE comconst_mod, ONLY: pi, dtvr |
---|
55 | |
---|
56 | PRIVATE |
---|
57 | PUBLIC :: etat0phys_netcdf |
---|
58 | |
---|
59 | include "iniprint.h" |
---|
60 | include "dimensions.h" |
---|
61 | include "paramet.h" |
---|
62 | include "comgeom2.h" |
---|
63 | include "dimsoil.h" |
---|
64 | include "clesphys.h" |
---|
65 | REAL, SAVE :: deg2rad |
---|
66 | REAL, SAVE, ALLOCATABLE :: tsol(:) |
---|
67 | INTEGER, SAVE :: iml_phys, jml_phys, llm_phys, ttm_phys, fid_phys |
---|
68 | REAL, ALLOCATABLE, SAVE :: lon_phys(:,:), lat_phys(:,:), levphys_ini(:) |
---|
69 | CHARACTER(LEN=256), PARAMETER :: oroparam="oro_params.nc" |
---|
70 | CHARACTER(LEN=256), PARAMETER :: orofname="Relief.nc", orogvar="RELIEF" |
---|
71 | CHARACTER(LEN=256), PARAMETER :: phyfname="ECPHY.nc", psrfvar="SP" |
---|
72 | CHARACTER(LEN=256), PARAMETER :: qsolvar="CDSW", tsrfvar="ST" |
---|
73 | |
---|
74 | |
---|
75 | CONTAINS |
---|
76 | |
---|
77 | |
---|
78 | !------------------------------------------------------------------------------- |
---|
79 | ! |
---|
80 | SUBROUTINE etat0phys_netcdf(masque, phis) |
---|
81 | ! |
---|
82 | !------------------------------------------------------------------------------- |
---|
83 | ! Purpose: Creates initial states |
---|
84 | !------------------------------------------------------------------------------- |
---|
85 | ! Notes: 1) This routine is designed to work for Earth |
---|
86 | ! 2) If masque(:,:)/=-99999., masque and phis are already known. |
---|
87 | ! Otherwise: compute it. |
---|
88 | !------------------------------------------------------------------------------- |
---|
89 | USE control_mod |
---|
90 | USE fonte_neige_mod |
---|
91 | USE pbl_surface_mod |
---|
92 | USE regr_horiz_time_climoz_m, ONLY: regr_horiz_time_climoz |
---|
93 | USE indice_sol_mod |
---|
94 | USE conf_phys_m, ONLY: conf_phys |
---|
95 | USE init_ssrf_m, ONLY: start_init_subsurf |
---|
96 | USE phys_state_var_mod, ONLY: beta_aridity, delta_tsurf, awake_dens, cv_gen, & |
---|
97 | ratqs_inter_ |
---|
98 | !use ioipsl_getincom |
---|
99 | IMPLICIT NONE |
---|
100 | !------------------------------------------------------------------------------- |
---|
101 | ! Arguments: |
---|
102 | REAL, INTENT(INOUT) :: masque(:,:) !--- Land mask dim(iip1,jjp1) |
---|
103 | REAL, INTENT(INOUT) :: phis (:,:) !--- Ground geopotential dim(iip1,jjp1) |
---|
104 | !------------------------------------------------------------------------------- |
---|
105 | ! Local variables: |
---|
106 | CHARACTER(LEN=256) :: modname="etat0phys_netcdf", fmt |
---|
107 | INTEGER :: i, j, l, ji, iml, jml |
---|
108 | LOGICAL :: read_mask |
---|
109 | REAL :: phystep, dummy |
---|
110 | REAL, DIMENSION(SIZE(masque,1),SIZE(masque,2)) :: masque_tmp,phiso |
---|
111 | REAL, DIMENSION(klon) :: sn, rugmer, run_off_lic_0, fder |
---|
112 | REAL, DIMENSION(klon,nbsrf) :: qsurf, snsrf |
---|
113 | REAL, DIMENSION(klon,nsoilmx,nbsrf) :: tsoil |
---|
114 | |
---|
115 | !--- Arguments for conf_phys |
---|
116 | LOGICAL :: ok_journe, ok_mensuel, ok_instan, ok_hf, ok_LES, callstats |
---|
117 | REAL :: solarlong0, seuil_inversion, fact_cldcon, facttemps |
---|
118 | LOGICAL :: ok_newmicro |
---|
119 | INTEGER :: iflag_radia, iflag_cldcon, iflag_ratqs |
---|
120 | REAL :: ratqsbas, ratqshaut, tau_ratqs |
---|
121 | LOGICAL :: ok_ade, ok_aie, ok_volcan, ok_alw, ok_cdnc, aerosol_couple, chemistry_couple |
---|
122 | INTEGER :: flag_aerosol |
---|
123 | INTEGER :: flag_aerosol_strat |
---|
124 | INTEGER :: flag_volc_surfstrat |
---|
125 | LOGICAL :: flag_aer_feedback |
---|
126 | LOGICAL :: flag_bc_internal_mixture |
---|
127 | REAL :: bl95_b0, bl95_b1 |
---|
128 | INTEGER :: read_climoz !--- Read ozone climatology |
---|
129 | REAL :: alp_offset |
---|
130 | LOGICAL :: filtre_oro=.false. |
---|
131 | |
---|
132 | INCLUDE "compbl.h" |
---|
133 | INCLUDE "alpale.h" |
---|
134 | |
---|
135 | deg2rad= pi/180.0 |
---|
136 | iml=assert_eq(SIZE(masque,1),SIZE(phis,1),TRIM(modname)//" iml") |
---|
137 | jml=assert_eq(SIZE(masque,2),SIZE(phis,2),TRIM(modname)//" jml") |
---|
138 | |
---|
139 | ! Physics configuration |
---|
140 | !******************************************************************************* |
---|
141 | CALL conf_phys( ok_journe, ok_mensuel, ok_instan, ok_hf, ok_LES, & |
---|
142 | callstats, & |
---|
143 | solarlong0,seuil_inversion, & |
---|
144 | fact_cldcon, facttemps,ok_newmicro,iflag_radia, & |
---|
145 | iflag_cldcon, & |
---|
146 | ratqsbas,ratqshaut,tau_ratqs, & |
---|
147 | ok_ade, ok_aie, ok_alw, ok_cdnc, ok_volcan, flag_volc_surfstrat, & |
---|
148 | aerosol_couple, chemistry_couple, flag_aerosol, flag_aerosol_strat, & |
---|
149 | flag_aer_feedback, flag_bc_internal_mixture, bl95_b0, bl95_b1, & |
---|
150 | read_climoz, alp_offset) |
---|
151 | CALL phys_state_var_init(read_climoz) |
---|
152 | |
---|
153 | !--- Initial atmospheric CO2 conc. from .def file |
---|
154 | co2_ppm0 = co2_ppm |
---|
155 | |
---|
156 | ! Compute ground geopotential, sub-cells quantities and possibly the mask. |
---|
157 | !******************************************************************************* |
---|
158 | read_mask=ANY(masque/=-99999.); masque_tmp=masque |
---|
159 | CALL start_init_orog(rlonv, rlatu, phis, masque_tmp) |
---|
160 | |
---|
161 | CALL getin('filtre_oro',filtre_oro) |
---|
162 | IF (filtre_oro) CALL filtreoro(size(phis,1),size(phis,2),phis,masque_tmp,rlatu) |
---|
163 | |
---|
164 | WRITE(fmt,"(i4,'i1)')")iml ; fmt='('//ADJUSTL(fmt) |
---|
165 | IF(.NOT.read_mask) THEN !--- Keep mask form orography |
---|
166 | masque=masque_tmp |
---|
167 | IF(prt_level>=1) THEN |
---|
168 | WRITE(lunout,*)'BUILT MASK :' |
---|
169 | WRITE(lunout,fmt) NINT(masque) |
---|
170 | END IF |
---|
171 | WHERE( masque(:,:)<EPSFRA) masque(:,:)=0. |
---|
172 | WHERE(1.-masque(:,:)<EPSFRA) masque(:,:)=1. |
---|
173 | END IF |
---|
174 | CALL gr_dyn_fi(1,iml,jml,klon,masque,zmasq) !--- Land mask to physical grid |
---|
175 | |
---|
176 | ! Compute tsol and qsol on physical grid, knowing phis on 2D grid. |
---|
177 | !******************************************************************************* |
---|
178 | CALL start_init_phys(rlonu, rlatv, phis) |
---|
179 | |
---|
180 | ! Some initializations. |
---|
181 | !******************************************************************************* |
---|
182 | sn (:) = 0.0 !--- Snow |
---|
183 | radsol(:) = 0.0 !--- Net radiation at ground |
---|
184 | rugmer(:) = 0.001 !--- Ocean rugosity |
---|
185 | !--- Ozone (zonal or 3D) interpolation in space and time (if 2nd arg is TRUE) |
---|
186 | IF(read_climoz>=1) CALL regr_horiz_time_climoz(read_climoz,ok_daily_climoz) |
---|
187 | |
---|
188 | ! Sub-surfaces initialization. |
---|
189 | !******************************************************************************* |
---|
190 | CALL start_init_subsurf(read_mask) |
---|
191 | |
---|
192 | ! Write physical initial state |
---|
193 | !******************************************************************************* |
---|
194 | WRITE(lunout,*)'phystep ',dtvr,iphysiq,nbapp_rad |
---|
195 | phystep = dtvr * FLOAT(iphysiq) |
---|
196 | radpas = NINT (86400./phystep/ FLOAT(nbapp_rad) ) |
---|
197 | WRITE(lunout,*)'phystep =', phystep, radpas |
---|
198 | |
---|
199 | ! Init: ftsol, snsrf, qsurf, tsoil, rain_fall, snow_fall, solsw, sollw, z0 |
---|
200 | !******************************************************************************* |
---|
201 | DO i=1,nbsrf; ftsol(:,i) = tsol; END DO |
---|
202 | DO i=1,nbsrf; snsrf(:,i) = sn; END DO |
---|
203 | falb_dir(:, :, is_ter) = 0.08 |
---|
204 | falb_dir(:, :, is_lic) = 0.6 |
---|
205 | falb_dir(:, :, is_oce) = 0.5 |
---|
206 | falb_dir(:, :, is_sic) = 0.6 |
---|
207 | |
---|
208 | !ym warning missing init for falb_dif => set to 0 |
---|
209 | falb_dif(:,:,:)=0 |
---|
210 | |
---|
211 | u10m(:,:)=0 |
---|
212 | v10m(:,:)=0 |
---|
213 | treedrg(:,:,:)=0 |
---|
214 | |
---|
215 | fevap(:,:) = 0. |
---|
216 | qsurf = 0. |
---|
217 | DO i=1,nbsrf; DO j=1,nsoilmx; tsoil(:,j,i) = tsol; END DO; END DO |
---|
218 | rain_fall = 0. |
---|
219 | snow_fall = 0. |
---|
220 | solsw = 165. |
---|
221 | solswfdiff = 1. |
---|
222 | sollw = -53. |
---|
223 | !ym warning missing init for sollwdown => set to 0 |
---|
224 | sollwdown = 0. |
---|
225 | t_ancien = 273.15 |
---|
226 | q_ancien = 0. |
---|
227 | ql_ancien = 0. |
---|
228 | qs_ancien = 0. |
---|
229 | prlw_ancien = 0. |
---|
230 | prsw_ancien = 0. |
---|
231 | prw_ancien = 0. |
---|
232 | agesno = 0. |
---|
233 | |
---|
234 | u_ancien = 0. |
---|
235 | v_ancien = 0. |
---|
236 | wake_delta_pbl_TKE(:,:,:)=0 |
---|
237 | wake_dens(:)=0 |
---|
238 | awake_dens = 0. |
---|
239 | cv_gen = 0. |
---|
240 | ale_bl = 0. |
---|
241 | ale_bl_trig =0. |
---|
242 | alp_bl=0. |
---|
243 | ale_wake=0. |
---|
244 | ale_bl_stat=0. |
---|
245 | |
---|
246 | cf_ancien = 0. |
---|
247 | rvc_ancien = 0. |
---|
248 | |
---|
249 | z0m(:,:)=0 ! ym missing 5th subsurface initialization |
---|
250 | |
---|
251 | z0m(:,is_oce) = rugmer(:) |
---|
252 | z0m(:,is_ter) = 0.01 !MAX(1.0e-05,zstd(:)*zsig(:)/2.0) |
---|
253 | z0m(:,is_lic) = 0.001 !MAX(1.0e-05,zstd(:)*zsig(:)/2.0) |
---|
254 | z0m(:,is_sic) = 0.001 |
---|
255 | z0h(:,:)=z0m(:,:) |
---|
256 | |
---|
257 | fder = 0.0 |
---|
258 | clwcon = 0.0 |
---|
259 | rnebcon = 0.0 |
---|
260 | ratqs = 0.0 |
---|
261 | run_off_lic_0 = 0.0 |
---|
262 | rugoro = 0.0 |
---|
263 | |
---|
264 | ! Before phyredem calling, surface modules and values to be saved in startphy.nc |
---|
265 | ! are initialized |
---|
266 | !******************************************************************************* |
---|
267 | dummy = 1.0 |
---|
268 | pbl_tke(:,:,:) = 1.e-8 |
---|
269 | zmax0(:) = 40. |
---|
270 | f0(:) = 1.e-5 |
---|
271 | sig1(:,:) = 0. |
---|
272 | w01(:,:) = 0. |
---|
273 | wake_deltat(:,:) = 0. |
---|
274 | wake_deltaq(:,:) = 0. |
---|
275 | wake_s(:) = 0. |
---|
276 | wake_cstar(:) = 0. |
---|
277 | wake_fip(:) = 0. |
---|
278 | wake_pe = 0. |
---|
279 | fm_therm = 0. |
---|
280 | entr_therm = 0. |
---|
281 | detr_therm = 0. |
---|
282 | awake_s = 0. |
---|
283 | |
---|
284 | CALL fonte_neige_init(run_off_lic_0) |
---|
285 | CALL pbl_surface_init( fder, snsrf, qsurf, tsoil ) |
---|
286 | |
---|
287 | IF (iflag_pbl>1 .AND. iflag_wake>=1 .AND. iflag_pbl_split >=1) then |
---|
288 | delta_tsurf = 0. |
---|
289 | beta_aridity = 0. |
---|
290 | end IF |
---|
291 | |
---|
292 | ratqs_inter_ = 0.002 |
---|
293 | CALL phyredem( "startphy.nc" ) |
---|
294 | |
---|
295 | ! WRITE(lunout,*)'CCCCCCCCCCCCCCCCCC REACTIVER SORTIE VISU DANS ETAT0' |
---|
296 | ! WRITE(lunout,*)'entree histclo' |
---|
297 | CALL histclo() |
---|
298 | |
---|
299 | END SUBROUTINE etat0phys_netcdf |
---|
300 | ! |
---|
301 | !------------------------------------------------------------------------------- |
---|
302 | |
---|
303 | |
---|
304 | !------------------------------------------------------------------------------- |
---|
305 | ! |
---|
306 | SUBROUTINE start_init_orog(lon_in,lat_in,phis,masque) |
---|
307 | ! |
---|
308 | !=============================================================================== |
---|
309 | ! Comment: |
---|
310 | ! This routine launch grid_noro, which computes parameters for SSO scheme as |
---|
311 | ! described in LOTT & MILLER (1997) and LOTT(1999). |
---|
312 | ! In case the file oroparam is present and the key read_orop is activated, |
---|
313 | ! grid_noro is bypassed and sub-cell parameters are read from the file. |
---|
314 | !=============================================================================== |
---|
315 | USE grid_noro_m, ONLY: grid_noro, read_noro |
---|
316 | USE logic_mod, ONLY: read_orop |
---|
317 | IMPLICIT NONE |
---|
318 | !------------------------------------------------------------------------------- |
---|
319 | ! Arguments: |
---|
320 | REAL, INTENT(IN) :: lon_in(:), lat_in(:) ! dim (iml) (jml) |
---|
321 | REAL, INTENT(INOUT) :: phis(:,:), masque(:,:) ! dim (iml,jml) |
---|
322 | !------------------------------------------------------------------------------- |
---|
323 | ! Local variables: |
---|
324 | CHARACTER(LEN=256) :: modname |
---|
325 | INTEGER :: fid, llm_tmp,ttm_tmp, iml,jml, iml_rel,jml_rel, itau(1) |
---|
326 | INTEGER :: ierr |
---|
327 | REAL :: lev(1), date, dt |
---|
328 | REAL, ALLOCATABLE :: lon_rad(:), lon_ini(:), lon_rel(:,:), relief_hi(:,:) |
---|
329 | REAL, ALLOCATABLE :: lat_rad(:), lat_ini(:), lat_rel(:,:), tmp_var (:,:) |
---|
330 | REAL, ALLOCATABLE :: zmea0(:,:), zstd0(:,:), zsig0(:,:) |
---|
331 | REAL, ALLOCATABLE :: zgam0(:,:), zthe0(:,:), zpic0(:,:), zval0(:,:) |
---|
332 | !------------------------------------------------------------------------------- |
---|
333 | modname="start_init_orog" |
---|
334 | iml=assert_eq(SIZE(lon_in),SIZE(phis,1),SIZE(masque,1),TRIM(modname)//" iml") |
---|
335 | jml=assert_eq(SIZE(lat_in),SIZE(phis,2),SIZE(masque,2),TRIM(modname)//" jml") |
---|
336 | |
---|
337 | !--- HIGH RESOLUTION OROGRAPHY |
---|
338 | CALL flininfo(orofname, iml_rel, jml_rel, llm_tmp, ttm_tmp, fid) |
---|
339 | |
---|
340 | ALLOCATE(lat_rel(iml_rel,jml_rel),lon_rel(iml_rel,jml_rel)) |
---|
341 | CALL flinopen(orofname, .FALSE., iml_rel, jml_rel, llm_tmp, lon_rel, lat_rel,& |
---|
342 | lev, ttm_tmp, itau, date, dt, fid) |
---|
343 | ALLOCATE(relief_hi(iml_rel,jml_rel)) |
---|
344 | CALL flinget(fid, orogvar, iml_rel, jml_rel, llm_tmp, ttm_tmp, 1,1, relief_hi) |
---|
345 | CALL flinclo(fid) |
---|
346 | |
---|
347 | !--- IF ANGLES ARE IN DEGREES, THEY ARE CONVERTED INTO RADIANS |
---|
348 | ALLOCATE(lon_ini(iml_rel),lat_ini(jml_rel)) |
---|
349 | lon_ini(:)=lon_rel(:,1); IF(MAXVAL(lon_rel)>pi) lon_ini=lon_ini*deg2rad |
---|
350 | lat_ini(:)=lat_rel(1,:); IF(MAXVAL(lat_rel)>pi) lat_ini=lat_ini*deg2rad |
---|
351 | |
---|
352 | !--- FIELDS ARE PROCESSED TO BE ON STANDARD ANGULAR DOMAINS |
---|
353 | ALLOCATE(lon_rad(iml_rel),lat_rad(jml_rel)) |
---|
354 | CALL conf_dat2d(orogvar, lon_ini, lat_ini, lon_rad, lat_rad, relief_hi,.FALSE.) |
---|
355 | DEALLOCATE(lon_ini,lat_ini) |
---|
356 | |
---|
357 | !--- COMPUTING THE REQUIRED FIELDS USING ROUTINE grid_noro |
---|
358 | WRITE(lunout,*) |
---|
359 | WRITE(lunout,*)'*** Compute parameters needed for gravity wave drag code ***' |
---|
360 | |
---|
361 | !--- ALLOCATIONS OF SUB-CELL SCALES QUANTITIES |
---|
362 | ALLOCATE(zmea0(iml,jml),zstd0(iml,jml)) !--- Mean orography and std deviation |
---|
363 | ALLOCATE(zsig0(iml,jml),zgam0(iml,jml)) !--- Slope and nisotropy |
---|
364 | zsig0(:,:)=0 !ym uninitialized variable |
---|
365 | zgam0(:,:)=0 !ym uninitialized variable |
---|
366 | ALLOCATE(zthe0(iml,jml)) !--- Highest slope orientation |
---|
367 | zthe0(:,:)=0 !ym uninitialized variable |
---|
368 | ALLOCATE(zpic0(iml,jml),zval0(iml,jml)) !--- Peaks and valley heights |
---|
369 | |
---|
370 | !--- READ SUB-CELL SCALES PARAMETERS FROM A FILE (AT RIGHT RESOLUTION) |
---|
371 | OPEN(UNIT=66,FILE=oroparam,STATUS='OLD',IOSTAT=ierr) |
---|
372 | IF(ierr==0.AND.read_orop) THEN |
---|
373 | CLOSE(UNIT=66) |
---|
374 | CALL read_noro(lon_in,lat_in,oroparam, & |
---|
375 | phis,zmea0,zstd0,zsig0,zgam0,zthe0,zpic0,zval0,masque) |
---|
376 | ELSE |
---|
377 | !--- CALL OROGRAPHY MODULE TO COMPUTE FIELDS |
---|
378 | CALL grid_noro(lon_rad,lat_rad,relief_hi,lon_in,lat_in, & |
---|
379 | phis,zmea0,zstd0,zsig0,zgam0,zthe0,zpic0,zval0,masque) |
---|
380 | END IF |
---|
381 | phis = phis * 9.81 |
---|
382 | phis(iml,:) = phis(1,:) |
---|
383 | DEALLOCATE(relief_hi,lon_rad,lat_rad) |
---|
384 | |
---|
385 | !--- PUT QUANTITIES TO PHYSICAL GRID |
---|
386 | CALL gr_dyn_fi(1,iml,jml,klon,zmea0,zmea); DEALLOCATE(zmea0) |
---|
387 | CALL gr_dyn_fi(1,iml,jml,klon,zstd0,zstd); DEALLOCATE(zstd0) |
---|
388 | CALL gr_dyn_fi(1,iml,jml,klon,zsig0,zsig); DEALLOCATE(zsig0) |
---|
389 | CALL gr_dyn_fi(1,iml,jml,klon,zgam0,zgam); DEALLOCATE(zgam0) |
---|
390 | CALL gr_dyn_fi(1,iml,jml,klon,zthe0,zthe); DEALLOCATE(zthe0) |
---|
391 | CALL gr_dyn_fi(1,iml,jml,klon,zpic0,zpic); DEALLOCATE(zpic0) |
---|
392 | CALL gr_dyn_fi(1,iml,jml,klon,zval0,zval); DEALLOCATE(zval0) |
---|
393 | |
---|
394 | |
---|
395 | END SUBROUTINE start_init_orog |
---|
396 | ! |
---|
397 | !------------------------------------------------------------------------------- |
---|
398 | |
---|
399 | |
---|
400 | !------------------------------------------------------------------------------- |
---|
401 | ! |
---|
402 | SUBROUTINE start_init_phys(lon_in,lat_in,phis) |
---|
403 | ! |
---|
404 | !=============================================================================== |
---|
405 | ! Purpose: Compute tsol and qsol, knowing phis. |
---|
406 | !=============================================================================== |
---|
407 | IMPLICIT NONE |
---|
408 | !------------------------------------------------------------------------------- |
---|
409 | ! Arguments: |
---|
410 | REAL, INTENT(IN) :: lon_in(:), lat_in(:) ! dim (iml) (jml2) |
---|
411 | REAL, INTENT(IN) :: phis(:,:) ! dim (iml,jml) |
---|
412 | !------------------------------------------------------------------------------- |
---|
413 | ! Local variables: |
---|
414 | CHARACTER(LEN=256) :: modname |
---|
415 | REAL :: date, dt |
---|
416 | INTEGER :: iml, jml, jml2, itau(1) |
---|
417 | REAL, ALLOCATABLE :: lon_rad(:), lon_ini(:), var_ana(:,:) |
---|
418 | REAL, ALLOCATABLE :: lat_rad(:), lat_ini(:) |
---|
419 | REAL, ALLOCATABLE :: ts(:,:), qs(:,:) |
---|
420 | !------------------------------------------------------------------------------- |
---|
421 | modname="start_init_phys" |
---|
422 | iml=assert_eq(SIZE(lon_in),SIZE(phis,1),TRIM(modname)//" iml") |
---|
423 | jml=SIZE(phis,2); jml2=SIZE(lat_in) |
---|
424 | |
---|
425 | WRITE(lunout,*)'Opening the surface analysis' |
---|
426 | CALL flininfo(phyfname, iml_phys, jml_phys, llm_phys, ttm_phys, fid_phys) |
---|
427 | WRITE(lunout,*) 'Values read: ', iml_phys, jml_phys, llm_phys, ttm_phys |
---|
428 | |
---|
429 | ALLOCATE(lat_phys(iml_phys,jml_phys),lon_phys(iml_phys,jml_phys)) |
---|
430 | ALLOCATE(levphys_ini(llm_phys)) |
---|
431 | CALL flinopen(phyfname, .FALSE., iml_phys, jml_phys, llm_phys, & |
---|
432 | lon_phys,lat_phys,levphys_ini,ttm_phys,itau,date,dt,fid_phys) |
---|
433 | |
---|
434 | !--- IF ANGLES ARE IN DEGREES, THEY ARE CONVERTED INTO RADIANS |
---|
435 | ALLOCATE(lon_ini(iml_phys),lat_ini(jml_phys)) |
---|
436 | lon_ini(:)=lon_phys(:,1); IF(MAXVAL(lon_phys)>pi) lon_ini=lon_ini*deg2rad |
---|
437 | lat_ini(:)=lat_phys(1,:); IF(MAXVAL(lat_phys)>pi) lat_ini=lat_ini*deg2rad |
---|
438 | |
---|
439 | ALLOCATE(var_ana(iml_phys,jml_phys),lon_rad(iml_phys),lat_rad(jml_phys)) |
---|
440 | CALL get_var_phys(tsrfvar,ts) !--- SURFACE TEMPERATURE |
---|
441 | CALL get_var_phys(qsolvar,qs) !--- SOIL MOISTURE |
---|
442 | CALL flinclo(fid_phys) |
---|
443 | DEALLOCATE(var_ana,lon_rad,lat_rad,lon_ini,lat_ini) |
---|
444 | |
---|
445 | !--- TSOL AND QSOL ON PHYSICAL GRID |
---|
446 | ALLOCATE(tsol(klon)) |
---|
447 | CALL gr_dyn_fi(1,iml,jml,klon,ts,tsol) |
---|
448 | CALL gr_dyn_fi(1,iml,jml,klon,qs,qsol) |
---|
449 | DEALLOCATE(ts,qs) |
---|
450 | |
---|
451 | CONTAINS |
---|
452 | |
---|
453 | !------------------------------------------------------------------------------- |
---|
454 | ! |
---|
455 | SUBROUTINE get_var_phys(title,field) |
---|
456 | ! |
---|
457 | !------------------------------------------------------------------------------- |
---|
458 | IMPLICIT NONE |
---|
459 | !------------------------------------------------------------------------------- |
---|
460 | ! Arguments: |
---|
461 | CHARACTER(LEN=*), INTENT(IN) :: title |
---|
462 | REAL, ALLOCATABLE, INTENT(INOUT) :: field(:,:) |
---|
463 | !------------------------------------------------------------------------------- |
---|
464 | ! Local variables: |
---|
465 | INTEGER :: tllm |
---|
466 | !------------------------------------------------------------------------------- |
---|
467 | SELECT CASE(title) |
---|
468 | CASE(psrfvar); tllm=0 |
---|
469 | CASE(tsrfvar,qsolvar); tllm=llm_phys |
---|
470 | END SELECT |
---|
471 | IF(ALLOCATED(field)) RETURN |
---|
472 | ALLOCATE(field(iml,jml)); field(:,:)=0. |
---|
473 | CALL flinget(fid_phys,title,iml_phys,jml_phys,tllm,ttm_phys,1,1,var_ana) |
---|
474 | CALL conf_dat2d(title, lon_ini, lat_ini, lon_rad, lat_rad, var_ana, .TRUE.) |
---|
475 | CALL interp_startvar(title, .TRUE., lon_rad, lat_rad, var_ana, & |
---|
476 | lon_in, lat_in, field) |
---|
477 | |
---|
478 | END SUBROUTINE get_var_phys |
---|
479 | ! |
---|
480 | !------------------------------------------------------------------------------- |
---|
481 | ! |
---|
482 | END SUBROUTINE start_init_phys |
---|
483 | ! |
---|
484 | !------------------------------------------------------------------------------- |
---|
485 | |
---|
486 | |
---|
487 | !------------------------------------------------------------------------------- |
---|
488 | ! |
---|
489 | SUBROUTINE interp_startvar(nam,ibeg,lon,lat,vari,lon2,lat2,varo) |
---|
490 | ! |
---|
491 | !------------------------------------------------------------------------------- |
---|
492 | USE inter_barxy_m, ONLY: inter_barxy |
---|
493 | IMPLICIT NONE |
---|
494 | !------------------------------------------------------------------------------- |
---|
495 | ! Arguments: |
---|
496 | CHARACTER(LEN=*), INTENT(IN) :: nam |
---|
497 | LOGICAL, INTENT(IN) :: ibeg |
---|
498 | REAL, INTENT(IN) :: lon(:), lat(:) ! dim (ii) (jj) |
---|
499 | REAL, INTENT(IN) :: vari(:,:) ! dim (ii,jj) |
---|
500 | REAL, INTENT(IN) :: lon2(:), lat2(:) ! dim (i1) (j2) |
---|
501 | REAL, INTENT(OUT) :: varo(:,:) ! dim (i1) (j1) |
---|
502 | !------------------------------------------------------------------------------- |
---|
503 | ! Local variables: |
---|
504 | CHARACTER(LEN=256) :: modname |
---|
505 | INTEGER :: ii, jj, i1, j1, j2 |
---|
506 | REAL, ALLOCATABLE :: vtmp(:,:) |
---|
507 | !------------------------------------------------------------------------------- |
---|
508 | modname="interp_startvar" |
---|
509 | ii=assert_eq(SIZE(lon), SIZE(vari,1),TRIM(modname)//" ii") |
---|
510 | jj=assert_eq(SIZE(lat), SIZE(vari,2),TRIM(modname)//" jj") |
---|
511 | i1=assert_eq(SIZE(lon2),SIZE(varo,1),TRIM(modname)//" i1") |
---|
512 | j1=SIZE(varo,2); j2=SIZE(lat2) |
---|
513 | ALLOCATE(vtmp(i1-1,j1)) |
---|
514 | IF(ibeg.AND.prt_level>1) THEN |
---|
515 | WRITE(lunout,*)"--------------------------------------------------------" |
---|
516 | WRITE(lunout,*)"$$$ Interpolation barycentrique pour "//TRIM(nam)//" $$$" |
---|
517 | WRITE(lunout,*)"--------------------------------------------------------" |
---|
518 | END IF |
---|
519 | CALL inter_barxy(lon, lat(:jj-1), vari, lon2(:i1-1), lat2, vtmp) |
---|
520 | CALL gr_int_dyn(vtmp, varo, i1-1, j1) |
---|
521 | |
---|
522 | END SUBROUTINE interp_startvar |
---|
523 | ! |
---|
524 | !------------------------------------------------------------------------------- |
---|
525 | ! |
---|
526 | !******************************************************************************* |
---|
527 | |
---|
528 | SUBROUTINE filtreoro(imp1,jmp1,phis,masque,rlatu) |
---|
529 | |
---|
530 | IMPLICIT NONE |
---|
531 | |
---|
532 | INTEGER imp1,jmp1 |
---|
533 | REAL, DIMENSION(imp1,jmp1) :: phis,masque |
---|
534 | REAL, DIMENSION(jmp1) :: rlatu |
---|
535 | REAL, DIMENSION(imp1) :: wwf |
---|
536 | REAL, DIMENSION(imp1,jmp1) :: phiso |
---|
537 | INTEGER :: ifiltre,ifi,ii,i,j |
---|
538 | REAL :: coslat0,ssz |
---|
539 | |
---|
540 | coslat0=0.5 |
---|
541 | phiso=phis |
---|
542 | do j=2,jmp1-1 |
---|
543 | print*,'avant if ',cos(rlatu(j)),coslat0 |
---|
544 | if (cos(rlatu(j))<coslat0) then |
---|
545 | ! nb de pts affectes par le filtrage de part et d'autre du pt |
---|
546 | ifiltre=(coslat0/cos(rlatu(j))-1.)/2. |
---|
547 | wwf=0. |
---|
548 | do i=1,ifiltre |
---|
549 | wwf(i)=1. |
---|
550 | enddo |
---|
551 | wwf(ifiltre+1)=(coslat0/cos(rlatu(j))-1.)/2.-ifiltre |
---|
552 | do i=1,imp1-1 |
---|
553 | if (masque(i,j)>0.9) then |
---|
554 | ssz=phis(i,j) |
---|
555 | do ifi=1,ifiltre+1 |
---|
556 | ii=i+ifi |
---|
557 | if (ii>imp1-1) ii=ii-imp1+1 |
---|
558 | ssz=ssz+wwf(ifi)*phis(ii,j) |
---|
559 | ii=i-ifi |
---|
560 | if (ii<1) ii=ii+imp1-1 |
---|
561 | ssz=ssz+wwf(ifi)*phis(ii,j) |
---|
562 | enddo |
---|
563 | phis(i,j)=ssz*cos(rlatu(j))/coslat0 |
---|
564 | endif |
---|
565 | enddo |
---|
566 | print*,'j=',j,coslat0/cos(rlatu(j)), (1.+2.*sum(wwf))*cos(rlatu(j))/coslat0 |
---|
567 | endif |
---|
568 | enddo |
---|
569 | call dump2d(imp1,jmp1,phis,'phis ') |
---|
570 | call dump2d(imp1,jmp1,masque,'masque ') |
---|
571 | call dump2d(imp1,jmp1,phis-phiso,'dphis ') |
---|
572 | |
---|
573 | END SUBROUTINE filtreoro |
---|
574 | |
---|
575 | |
---|
576 | END MODULE etat0phys |
---|