source: LMDZ6/branches/Amaury_dev/libf/phylmd/dyn1d/mod_1D_amma_read.F90

Last change on this file was 5158, checked in by abarral, 5 months ago

Add missing klon on strataer_emiss_mod.F90
Correct various missing explicit declarations
Replace tabs by spaces (tabs are not part of the fortran charset)
Continue cleaning modules
Removed unused arguments and variables
File size: 13.5 KB
Line 
1MODULE mod_1D_amma_read
2  USE netcdf, ONLY: nf90_get_var, nf90_open, nf90_noerr, nf90_open, nf90_nowrite, &
3          nf90_inq_dimid, nf90_inquire_dimension, nf90_strerror, nf90_inq_varid
4  !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
5  !Declarations specifiques au cas AMMA
6  CHARACTER*80 :: fich_amma
7  ! Option du cas AMMA ou on impose la discretisation verticale (Ap,Bp)
8  INTEGER nlev_amma, nt_amma
9
10  INTEGER year_ini_amma, day_ini_amma, mth_ini_amma
11  REAL heure_ini_amma
12  REAL day_ju_ini_amma   ! Julian day of amma first day
13  parameter (year_ini_amma = 2006)
14  parameter (mth_ini_amma = 7)
15  parameter (day_ini_amma = 10)  ! 10 = 10Juil2006
16  parameter (heure_ini_amma = 0.) !0h en secondes
17  REAL dt_amma
18  parameter (dt_amma = 1800.)
19
20  !profils initiaux:
21  REAL, ALLOCATABLE :: plev_amma(:)
22
23  REAL, ALLOCATABLE :: z_amma(:)
24  REAL, ALLOCATABLE :: th_amma(:), q_amma(:)
25  REAL, ALLOCATABLE :: u_amma(:)
26  REAL, ALLOCATABLE :: v_amma(:)
27
28  REAL, ALLOCATABLE :: th_ammai(:), q_ammai(:)
29  REAL, ALLOCATABLE :: u_ammai(:)
30  REAL, ALLOCATABLE :: v_ammai(:)
31  REAL, ALLOCATABLE :: vitw_ammai(:)
32  REAL, ALLOCATABLE :: ht_ammai(:)
33  REAL, ALLOCATABLE :: hq_ammai(:)
34  REAL, ALLOCATABLE :: vt_ammai(:)
35  REAL, ALLOCATABLE :: vq_ammai(:)
36
37  !forcings
38  REAL, ALLOCATABLE :: ht_amma(:, :)
39  REAL, ALLOCATABLE :: hq_amma(:, :)
40  REAL, ALLOCATABLE :: vitw_amma(:, :)
41  REAL, ALLOCATABLE :: lat_amma(:), sens_amma(:)
42
43  !champs interpoles
44  REAL, ALLOCATABLE :: vitw_profamma(:)
45  REAL, ALLOCATABLE :: ht_profamma(:)
46  REAL, ALLOCATABLE :: hq_profamma(:)
47  REAL lat_profamma, sens_profamma
48  REAL, ALLOCATABLE :: vt_profamma(:)
49  REAL, ALLOCATABLE :: vq_profamma(:)
50  REAL, ALLOCATABLE :: th_profamma(:)
51  REAL, ALLOCATABLE :: q_profamma(:)
52  REAL, ALLOCATABLE :: u_profamma(:)
53  REAL, ALLOCATABLE :: v_profamma(:)
54
55
56CONTAINS
57
58  SUBROUTINE read_1D_cases
59    IMPLICIT NONE
60
61    INTEGER nid, rid, ierr
62
63    fich_amma = 'amma.nc'
64    PRINT*, 'fich_amma ', fich_amma
65    ierr = nf90_open(fich_amma, nf90_nowrite, nid)
66    PRINT*, 'fich_amma,nf90_nowrite,nid ', fich_amma, nf90_nowrite, nid
67    IF (ierr/=nf90_noerr) THEN
68      WRITE(*, *) 'ERROR: GROS Pb opening forcings nc file '
69      WRITE(*, *) nf90_strerror(ierr)
70      stop ""
71    endif
72    !.......................................................................
73    ierr = nf90_inq_dimid(nid, 'lev', rid)
74    IF (ierr/=nf90_noerr) THEN
75      PRINT*, 'Oh probleme lecture dimension zz'
76    ENDIF
77    ierr = nf90_inquire_dimension(nid, rid, len = nlev_amma)
78    PRINT*, 'OK nid,rid,nlev_amma', nid, rid, nlev_amma
79    !.......................................................................
80    ierr = nf90_inq_dimid(nid, 'time', rid)
81    PRINT*, 'nid,rid', nid, rid
82    nt_amma = 0
83    IF (ierr/=nf90_noerr) THEN
84      stop 'probleme lecture dimension sens'
85    ENDIF
86    ierr = nf90_inquire_dimension(nid, rid, len = nt_amma)
87    PRINT*, 'nid,rid,nlev_amma', nid, rid, nt_amma
88
89    !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
90    !profils initiaux:
91    allocate(plev_amma(nlev_amma))
92
93    allocate(z_amma(nlev_amma))
94    allocate(th_amma(nlev_amma), q_amma(nlev_amma))
95    allocate(u_amma(nlev_amma))
96    allocate(v_amma(nlev_amma))
97
98    !forcings
99    allocate(ht_amma(nlev_amma, nt_amma))
100    allocate(hq_amma(nlev_amma, nt_amma))
101    allocate(vitw_amma(nlev_amma, nt_amma))
102    allocate(lat_amma(nt_amma), sens_amma(nt_amma))
103
104    !profils initiaux:
105    allocate(th_ammai(nlev_amma), q_ammai(nlev_amma))
106    allocate(u_ammai(nlev_amma))
107    allocate(v_ammai(nlev_amma))
108    allocate(vitw_ammai(nlev_amma))
109    allocate(ht_ammai(nlev_amma))
110    allocate(hq_ammai(nlev_amma))
111    allocate(vt_ammai(nlev_amma))
112    allocate(vq_ammai(nlev_amma))
113
114    !champs interpoles
115    allocate(vitw_profamma(nlev_amma))
116    allocate(ht_profamma(nlev_amma))
117    allocate(hq_profamma(nlev_amma))
118    allocate(vt_profamma(nlev_amma))
119    allocate(vq_profamma(nlev_amma))
120    allocate(th_profamma(nlev_amma))
121    allocate(q_profamma(nlev_amma))
122    allocate(u_profamma(nlev_amma))
123    allocate(v_profamma(nlev_amma))
124
125    PRINT*, 'Allocations OK'
126    CALL read_amma(nid, nlev_amma, nt_amma                                  &
127            , z_amma, plev_amma, th_amma, q_amma, u_amma, v_amma, vitw_amma         &
128            , ht_amma, hq_amma, sens_amma, lat_amma)
129
130  END SUBROUTINE read_1D_cases
131
132
133  !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
134  SUBROUTINE deallocate_1D_cases
135    !profils initiaux:
136    deallocate(plev_amma)
137
138    deallocate(z_amma)
139    deallocate(th_amma, q_amma)
140    deallocate(u_amma)
141    deallocate(v_amma)
142
143    deallocate(th_ammai, q_ammai)
144    deallocate(u_ammai)
145    deallocate(v_ammai)
146    deallocate(vitw_ammai)
147    deallocate(ht_ammai)
148    deallocate(hq_ammai)
149    deallocate(vt_ammai)
150    deallocate(vq_ammai)
151
152    !forcings
153    deallocate(ht_amma)
154    deallocate(hq_amma)
155    deallocate(vitw_amma)
156    deallocate(lat_amma, sens_amma)
157
158    !champs interpoles
159    deallocate(vitw_profamma)
160    deallocate(ht_profamma)
161    deallocate(hq_profamma)
162    deallocate(vt_profamma)
163    deallocate(vq_profamma)
164    deallocate(th_profamma)
165    deallocate(q_profamma)
166    deallocate(u_profamma)
167    deallocate(v_profamma)
168  END SUBROUTINE deallocate_1D_cases
169
170
171  !=====================================================================
172  SUBROUTINE read_amma(nid, nlevel, ntime                          &
173          , zz, pp, temp, qv, u, v, dw                   &
174          , dt, dq, sens, flat)
175
176    !program reading forcings of the AMMA case study
177    IMPLICIT NONE
178
179    INTEGER ntime, nlevel
180
181    REAL zz(nlevel)
182    REAL temp(nlevel), pp(nlevel)
183    REAL qv(nlevel), u(nlevel)
184    REAL v(nlevel)
185    REAL dw(nlevel, ntime)
186    REAL dt(nlevel, ntime)
187    REAL dq(nlevel, ntime)
188    REAL flat(ntime), sens(ntime)
189
190    INTEGER nid, ierr, rid
191    INTEGER nbvar3d
192    parameter(nbvar3d = 30)
193    INTEGER var3didin(nbvar3d)
194
195    ierr = nf90_inq_varid(nid, "zz", var3didin(1))
196    IF(ierr/=nf90_noerr) THEN
197      WRITE(*, *) nf90_strerror(ierr)
198      stop 'lev'
199    endif
200
201    ierr = nf90_inq_varid(nid, "temp", var3didin(2))
202    IF(ierr/=nf90_noerr) THEN
203      WRITE(*, *) nf90_strerror(ierr)
204      stop 'temp'
205    endif
206
207    ierr = nf90_inq_varid(nid, "qv", var3didin(3))
208    IF(ierr/=nf90_noerr) THEN
209      WRITE(*, *) nf90_strerror(ierr)
210      stop 'qv'
211    endif
212
213    ierr = nf90_inq_varid(nid, "u", var3didin(4))
214    IF(ierr/=nf90_noerr) THEN
215      WRITE(*, *) nf90_strerror(ierr)
216      stop 'u'
217    endif
218
219    ierr = nf90_inq_varid(nid, "v", var3didin(5))
220    IF(ierr/=nf90_noerr) THEN
221      WRITE(*, *) nf90_strerror(ierr)
222      stop 'v'
223    endif
224
225    ierr = nf90_inq_varid(nid, "dw", var3didin(6))
226    IF(ierr/=nf90_noerr) THEN
227      WRITE(*, *) nf90_strerror(ierr)
228      stop 'dw'
229    endif
230
231    ierr = nf90_inq_varid(nid, "dt", var3didin(7))
232    IF(ierr/=nf90_noerr) THEN
233      WRITE(*, *) nf90_strerror(ierr)
234      stop 'dt'
235    endif
236
237    ierr = nf90_inq_varid(nid, "dq", var3didin(8))
238    IF(ierr/=nf90_noerr) THEN
239      WRITE(*, *) nf90_strerror(ierr)
240      stop 'dq'
241    endif
242
243    ierr = nf90_inq_varid(nid, "sens", var3didin(9))
244    IF(ierr/=nf90_noerr) THEN
245      WRITE(*, *) nf90_strerror(ierr)
246      stop 'sens'
247    endif
248
249    ierr = nf90_inq_varid(nid, "flat", var3didin(10))
250    IF(ierr/=nf90_noerr) THEN
251      WRITE(*, *) nf90_strerror(ierr)
252      stop 'flat'
253    endif
254
255    ierr = nf90_inq_varid(nid, "pp", var3didin(11))
256    IF(ierr/=nf90_noerr) THEN
257      WRITE(*, *) nf90_strerror(ierr)
258    endif
259
260    !dimensions lecture
261    !      CALL catchaxis(nid,ntime,nlevel,time,z,ierr)
262
263    ierr = nf90_get_var(nid, var3didin(1), zz)
264    IF(ierr/=nf90_noerr) THEN
265      WRITE(*, *) nf90_strerror(ierr)
266      stop "getvarup"
267    endif
268    !          WRITE(*,*)'lecture z ok',zz
269
270    ierr = nf90_get_var(nid, var3didin(2), temp)
271    IF(ierr/=nf90_noerr) THEN
272      WRITE(*, *) nf90_strerror(ierr)
273      stop "getvarup"
274    endif
275    !          WRITE(*,*)'lecture th ok',temp
276
277    ierr = nf90_get_var(nid, var3didin(3), qv)
278    IF(ierr/=nf90_noerr) THEN
279      WRITE(*, *) nf90_strerror(ierr)
280      stop "getvarup"
281    endif
282    !          WRITE(*,*)'lecture qv ok',qv
283
284    ierr = nf90_get_var(nid, var3didin(4), u)
285    IF(ierr/=nf90_noerr) THEN
286      WRITE(*, *) nf90_strerror(ierr)
287      stop "getvarup"
288    endif
289    !          WRITE(*,*)'lecture u ok',u
290
291    ierr = nf90_get_var(nid, var3didin(5), v)
292    IF(ierr/=nf90_noerr) THEN
293      WRITE(*, *) nf90_strerror(ierr)
294      stop "getvarup"
295    endif
296    !          WRITE(*,*)'lecture v ok',v
297
298    ierr = nf90_get_var(nid, var3didin(6), dw)
299    IF(ierr/=nf90_noerr) THEN
300      WRITE(*, *) nf90_strerror(ierr)
301      stop "getvarup"
302    endif
303    !          WRITE(*,*)'lecture w ok',dw
304
305    ierr = nf90_get_var(nid, var3didin(7), dt)
306    IF(ierr/=nf90_noerr) THEN
307      WRITE(*, *) nf90_strerror(ierr)
308      stop "getvarup"
309    endif
310    !          WRITE(*,*)'lecture dt ok',dt
311
312    ierr = nf90_get_var(nid, var3didin(8), dq)
313    IF(ierr/=nf90_noerr) THEN
314      WRITE(*, *) nf90_strerror(ierr)
315      stop "getvarup"
316    endif
317    !          WRITE(*,*)'lecture dq ok',dq
318
319    ierr = nf90_get_var(nid, var3didin(9), sens)
320    IF(ierr/=nf90_noerr) THEN
321      WRITE(*, *) nf90_strerror(ierr)
322      stop "getvarup"
323    endif
324    !          WRITE(*,*)'lecture sens ok',sens
325
326    ierr = nf90_get_var(nid, var3didin(10), flat)
327    IF(ierr/=nf90_noerr) THEN
328      WRITE(*, *) nf90_strerror(ierr)
329      stop "getvarup"
330    endif
331    !          WRITE(*,*)'lecture flat ok',flat
332
333    ierr = nf90_get_var(nid, var3didin(11), pp)
334    IF(ierr/=nf90_noerr) THEN
335      WRITE(*, *) nf90_strerror(ierr)
336      stop "getvarup"
337    endif
338    !          WRITE(*,*)'lecture pp ok',pp
339
340  END SUBROUTINE  read_amma
341  !======================================================================
342  SUBROUTINE interp_amma_time(day, day1, annee_ref                     &
343          , year_ini_amma, day_ini_amma, nt_amma, dt_amma, nlev_amma       &
344          , vitw_amma, ht_amma, hq_amma, lat_amma, sens_amma               &
345          , vitw_prof, ht_prof, hq_prof, lat_prof, sens_prof)
346
347    USE lmdz_compar1d
348
349    IMPLICIT NONE
350
351    !---------------------------------------------------------------------------------------
352    ! Time interpolation of a 2D field to the timestep corresponding to day
353
354    ! day: current julian day (e.g. 717538.2)
355    ! day1: first day of the simulation
356    ! nt_amma: total nb of data in the forcing (e.g. 48 for AMMA)
357    ! dt_amma: total time interval (in sec) between 2 forcing data (e.g. 30min for AMMA)
358    !---------------------------------------------------------------------------------------
359
360    ! inputs:
361    INTEGER annee_ref
362    INTEGER nt_amma, nlev_amma
363    INTEGER year_ini_amma
364    REAL day, day1, day_ini_amma, dt_amma
365    REAL vitw_amma(nlev_amma, nt_amma)
366    REAL ht_amma(nlev_amma, nt_amma)
367    REAL hq_amma(nlev_amma, nt_amma)
368    REAL lat_amma(nt_amma)
369    REAL sens_amma(nt_amma)
370    ! outputs:
371    REAL vitw_prof(nlev_amma)
372    REAL ht_prof(nlev_amma)
373    REAL hq_prof(nlev_amma)
374    REAL lat_prof, sens_prof
375    ! local:
376    INTEGER it_amma1, it_amma2, k
377    REAL timeit, time_amma1, time_amma2, frac
378
379    IF (forcing_type==6) THEN
380      ! Check that initial day of the simulation consistent with AMMA case:
381      IF (annee_ref/=2006) THEN
382        PRINT*, 'Pour AMMA, annee_ref doit etre 2006'
383        PRINT*, 'Changer annee_ref dans run.def'
384        stop
385      endif
386      IF (annee_ref==2006 .AND. day1<day_ini_amma) THEN
387        PRINT*, 'AMMA a débuté le 10 juillet 2006', day1, day_ini_amma
388        PRINT*, 'Changer dayref dans run.def'
389        stop
390      endif
391      IF (annee_ref==2006 .AND. day1>day_ini_amma + 1) THEN
392        PRINT*, 'AMMA a fini le 11 juillet'
393        PRINT*, 'Changer dayref ou nday dans run.def'
394        stop
395      endif
396    endif
397
398    ! Determine timestep relative to the 1st day of AMMA:
399    !       timeit=(day-day1)*86400.
400    !       if (annee_ref.EQ.1992) THEN
401    !        timeit=(day-day_ini_toga)*86400.
402    !       else
403    !        timeit=(day+61.-1.)*86400. ! 61 days between Nov01 and Dec31 1992
404    !       endif
405    timeit = (day - day_ini_amma) * 86400
406
407    ! Determine the closest observation times:
408    !       it_amma1=INT(timeit/dt_amma)+1
409    !       it_amma2=it_amma1 + 1
410    !       time_amma1=(it_amma1-1)*dt_amma
411    !       time_amma2=(it_amma2-1)*dt_amma
412
413    it_amma1 = INT(timeit / dt_amma) + 1
414    IF (it_amma1 == nt_amma) THEN
415      it_amma2 = it_amma1
416    ELSE
417      it_amma2 = it_amma1 + 1
418    ENDIF
419    time_amma1 = (it_amma1 - 1) * dt_amma
420    time_amma2 = (it_amma2 - 1) * dt_amma
421
422    IF (it_amma1 > nt_amma) THEN
423      WRITE(*, *) 'PB-stop: day, it_amma1, it_amma2, timeit: '            &
424              , day, day_ini_amma, it_amma1, it_amma2, timeit / 86400.
425      stop
426    endif
427
428    ! time interpolation:
429    IF (it_amma1 == it_amma2) THEN
430      frac = 0.
431    ELSE
432      frac = (time_amma2 - timeit) / (time_amma2 - time_amma1)
433      frac = max(frac, 0.0)
434    ENDIF
435
436    lat_prof = lat_amma(it_amma2)                                       &
437            - frac * (lat_amma(it_amma2) - lat_amma(it_amma1))
438    sens_prof = sens_amma(it_amma2)                                     &
439            - frac * (sens_amma(it_amma2) - sens_amma(it_amma1))
440
441    DO k = 1, nlev_amma
442      vitw_prof(k) = vitw_amma(k, it_amma2)                               &
443              - frac * (vitw_amma(k, it_amma2) - vitw_amma(k, it_amma1))
444      ht_prof(k) = ht_amma(k, it_amma2)                                   &
445              - frac * (ht_amma(k, it_amma2) - ht_amma(k, it_amma1))
446      hq_prof(k) = hq_amma(k, it_amma2)                                   &
447              - frac * (hq_amma(k, it_amma2) - hq_amma(k, it_amma1))
448    enddo
449
450    RETURN
451  END
452
453END MODULE mod_1D_amma_read
Note: See TracBrowser for help on using the repository browser.