GCC Code Coverage Report
Directory: ./ Exec Total Coverage
File: dynphy_lonlat/phylmd/etat0phys_netcdf.F90 Lines: 0 232 0.0 %
Date: 2023-06-30 12:51:15 Branches: 0 661 0.0 %

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