source: LMDZ5/trunk/libf/dyn3dmem/startvar.F @ 1632

Last change on this file since 1632 was 1632, checked in by Laurent Fairhead, 12 years ago

Import initial du répertoire dyn3dmem

Attention! ceci n'est qu'une version préliminaire du code "basse mémoire":
le code contenu dans ce répertoire est basé sur la r1320 et a donc besoin
d'être mis à jour par rapport à la dynamique parallèle d'aujourd'hui.
Ce code est toutefois mis à disposition pour circonvenir à des problèmes
de mémoire que certaines configurations du modèle pourraient rencontrer.
Dans l'état, il compile et tourne sur vargas et au CCRT

Initial import of dyn3dmem

Warning! this is just a preliminary version of the memory light code:
it is based on r1320 of the code and thus needs to be updated before
it can replace the present dyn3dpar code. It is nevertheless put at your
disposal to circumvent some memory problems some LMDZ configurations may
encounter. In its present state, it will compile and run on vargas and CCRT

File size: 40.4 KB
Line 
1!
2! $Id: startvar.F 1299 2010-01-20 14:27:21Z fairhead $
3!
4      MODULE startvar
5#ifdef CPP_EARTH
6! This module is designed to work for Earth (and with ioipsl)
7    !
8    !
9    !      There are three ways to access data from the database of atmospheric data which
10    !       can be used to initialize the model. This depends on the type of field which needs
11    !       to be extracted.
12    !       We will details the possible arguments to startget here :
13    !
14    !        - A 2D variable on the dynamical grid :
15    !           CALL startget_phys2d(varname, iml, jml, lon_in, lat_in, champ, val_ex, jml2, lon_in2, lat_in2, interbar )             
16    !
17    !        - A 1D variable on the physical grid :
18    !            CALL startget_phys1d(varname, iml, jml, lon_in, lat_in, nbindex, champ, val_exp, jml2, lon_in2, lat_in2, interbar )
19    !
20    !
21    !         - A 3D variable on the dynamical grid :
22    !            CALL startget_dyn(varname, iml, jml, lon_in, lat_in, lml, pls, workvar, champ, val_exp, jml2, lon_in2, lat_in2, interbar )
23    !
24    !
25    !         There is special constraint on the atmospheric data base except that the
26    !         the data needs to be in netCDF and the variables should have the the following
27    !        names in the file :
28    !
29    !      'RELIEF'  : High resolution orography
30    !       'ST'            : Surface temperature
31    !       'CDSW'     : Soil moisture
32    !       'Z'               : Surface geopotential
33    !       'SP'            : Surface pressure
34    !        'U'              : East ward wind
35    !        'V'              : Northward wind
36    !        'TEMP'             : Temperature
37    !        'R'             : Relative humidity
38    !     
39      USE ioipsl
40    !
41    !
42      IMPLICIT NONE
43    !
44    !
45      PRIVATE
46      public startget_phys2d, startget_phys1d, startget_dyn
47    !
48      INTEGER, SAVE :: fid_phys, fid_dyn
49      INTEGER, SAVE  :: iml_phys, iml_rel, iml_dyn
50      INTEGER, SAVE :: jml_phys,  jml_rel, jml_dyn
51      INTEGER, SAVE ::  llm_dyn, ttm_dyn
52      REAL, ALLOCATABLE, SAVE, DIMENSION (:,:)  :: lon_phys, lon_rug,
53     . lon_alb, lon_rel, lon_dyn
54      REAL, ALLOCATABLE, SAVE, DIMENSION (:,:)  :: lat_phys, lat_rug,
55     . lat_alb, lat_rel, lat_dyn
56      REAL, ALLOCATABLE, SAVE, DIMENSION (:)  :: levdyn_ini
57      REAL, ALLOCATABLE, SAVE, DIMENSION (:,:)  :: relief, zstd, zsig,
58     . zgam, zthe, zpic, zval
59      REAL, ALLOCATABLE, SAVE, DIMENSION (:,:)  :: rugo, masque, phis
60    !
61      REAL, ALLOCATABLE, SAVE, DIMENSION (:,:)  :: tsol, qsol, psol_dyn
62      REAL, ALLOCATABLE, SAVE, DIMENSION (:,:,:)  ::   var_ana3d
63    !
64      CONTAINS
65    !
66    !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
67    !
68      SUBROUTINE startget_phys2d(varname, iml, jml, lon_in, lat_in,
69     . champ, val_exp, jml2, lon_in2, lat_in2 , interbar, masque_lu )
70    !
71    !    There is a big mess with the size in logitude, should it be iml or iml+1.
72    !    I have chosen to use the iml+1 as an argument to this routine and we declare
73    !   internaly smaler fields when needed. This needs to be cleared once and for all in LMDZ.
74    !  A convention is required.
75    !
76    !
77      CHARACTER*(*), INTENT(in) :: varname
78      INTEGER, INTENT(in) :: iml, jml ,jml2
79      REAL, INTENT(in) :: lon_in(iml), lat_in(jml)
80      REAL, INTENT(in) :: lon_in2(iml), lat_in2(jml2)
81      REAL, INTENT(inout) :: champ(iml,jml)
82      REAL, INTENT(in) :: val_exp
83      REAL, INTENT(in), optional :: masque_lu(iml,jml)
84      LOGICAL interbar
85    !
86    !   This routine only works if the variable does not exist or is constant
87    !
88      IF ( MINVAL(champ(:,:)).EQ.MAXVAL(champ(:,:)) .AND.
89     .MINVAL(champ(:,:)).EQ.val_exp ) THEN
90          !
91          SELECTCASE(varname)
92              !
93              CASE ('relief')
94                  !
95                  !  If we do not have the orography we need to get it
96                  !
97                  IF ( .NOT.ALLOCATED(relief)) THEN
98                      !
99                    if (present(masque_lu)) then
100                      CALL start_init_orog( iml, jml, lon_in, lat_in,
101     .                    jml2,lon_in2,lat_in2, interbar, masque_lu )
102                    else
103                      CALL start_init_orog( iml, jml, lon_in, lat_in,
104     .                    jml2,lon_in2,lat_in2, interbar)
105                    endif
106                      !
107                  ENDIF
108                  !
109                  IF (SIZE(relief) .NE. SIZE(champ)) THEN
110                      !
111                      WRITE(*,*) 'STARTVAR module has been',
112     .' initialized to the wrong size'
113                      STOP
114                      !
115                  ENDIF
116                  !
117                  champ(:,:) = relief(:,:)
118                  !
119              CASE ('rugosite')
120                  !
121                  !  If we do not have the orography we need to get it
122                  !
123                  IF ( .NOT.ALLOCATED(rugo)) THEN
124                      !
125                      CALL start_init_orog( iml, jml, lon_in, lat_in,
126     .                    jml2,lon_in2,lat_in2 , interbar )
127                      !
128                  ENDIF
129                  !
130                  IF (SIZE(rugo) .NE. SIZE(champ)) THEN
131                      !
132                      WRITE(*,*)
133     .  'STARTVAR module has been initialized to the wrong size'
134                      STOP
135                      !
136                  ENDIF
137                  !
138                  champ(:,:) = rugo(:,:)
139                  !
140              CASE ('masque')
141                  !
142                  !  If we do not have the orography we need to get it
143                  !
144                  IF ( .NOT.ALLOCATED(masque)) THEN
145                      !
146                      CALL start_init_orog( iml, jml, lon_in, lat_in,
147     .                     jml2,lon_in2,lat_in2 , interbar )
148                      !
149                  ENDIF
150                  !
151                  IF (SIZE(masque) .NE. SIZE(champ)) THEN
152                      !
153                      WRITE(*,*)
154     .   'STARTVAR module has been initialized to the wrong size'
155                      STOP
156                      !
157                  ENDIF
158                  !
159                  champ(:,:) = masque(:,:)
160                  !
161              CASE ('surfgeo')
162                  !
163                  !  If we do not have the orography we need to get it
164                  !
165                  IF ( .NOT.ALLOCATED(phis)) THEN
166                      !
167                      CALL start_init_orog( iml, jml, lon_in, lat_in,
168     .                   jml2,lon_in2, lat_in2 , interbar )
169                      !
170                  ENDIF
171                  !
172                  IF (SIZE(phis) .NE. SIZE(champ)) THEN
173                      !
174                      WRITE(*,*)
175     . 'STARTVAR module has been initialized to the wrong size'
176                      STOP
177                      !
178                  ENDIF
179                  !
180                  champ(:,:) = phis(:,:)
181                  !
182              CASE ('psol')
183                  !
184                  !  If we do not have the orography we need to get it
185                  !
186                  IF ( .NOT.ALLOCATED(psol_dyn)) THEN
187                      !
188                      CALL start_init_dyn( iml, jml, lon_in, lat_in,
189     .                   jml2,lon_in2, lat_in2 , interbar )
190                      !
191                  ENDIF
192                  !
193                  IF (SIZE(psol_dyn) .NE. SIZE(champ)) THEN
194                      !
195                      WRITE(*,*)
196     . 'STARTVAR module has been initialized to the wrong size'
197                      STOP
198                      !
199                  ENDIF
200                  !
201                  champ(:,:) = psol_dyn(:,:)
202                  !
203              CASE DEFAULT
204                  !
205                  WRITE(*,*) 'startget_phys2d'
206                  WRITE(*,*) 'No rule is present to extract variable',
207     .                 varname(:LEN_TRIM(varname)),' from any data set'
208                  STOP
209                  !
210          END SELECT
211          !
212      ELSE
213          !
214          ! There are a few fields we might need if we need to interpolate 3D filed. Thus if they come through here we
215          ! will catch them
216          !
217          SELECTCASE(varname)
218              !
219              CASE ('surfgeo')
220                  !
221                  IF ( .NOT.ALLOCATED(phis)) THEN
222                      ALLOCATE(phis(iml,jml))
223                  ENDIF
224                  !
225                  IF (SIZE(phis) .NE. SIZE(champ)) THEN
226                      !
227                      WRITE(*,*)
228     .  'STARTVAR module has been initialized to the wrong size'
229                      STOP
230                      !
231                  ENDIF
232                  !
233                  phis(:,:) = champ(:,:)
234                  !
235          END SELECT
236          !
237      ENDIF
238    !
239      END SUBROUTINE startget_phys2d
240    !
241    !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
242    !
243      SUBROUTINE start_init_orog ( iml,jml,lon_in, lat_in,jml2,lon_in2 ,
244     ,   lat_in2 , interbar, masque_lu )
245    !
246      INTEGER, INTENT(in) :: iml, jml, jml2
247      REAL, INTENT(in) :: lon_in(iml), lat_in(jml)
248      REAL, INTENT(in) :: lon_in2(iml), lat_in2(jml2)
249      REAL, intent(in), optional :: masque_lu(iml,jml)
250      LOGICAL interbar
251    !
252    !  LOCAL
253    !
254      LOGICAL interbar2
255      REAL :: lev(1), date, dt,chmin,chmax
256      INTEGER :: itau(1), fid
257      INTEGER ::  llm_tmp, ttm_tmp
258      INTEGER :: i, j
259      INTEGER :: iret
260      CHARACTER*25 title
261      REAL, ALLOCATABLE :: relief_hi(:,:)
262      REAL, ALLOCATABLE :: lon_rad(:), lat_rad(:)
263      REAL, ALLOCATABLE :: lon_ini(:), lat_ini(:)
264      REAL, ALLOCATABLE :: tmp_var(:,:)
265      INTEGER, ALLOCATABLE :: tmp_int(:,:)
266    !
267      CHARACTER*120 :: orogfname
268      LOGICAL :: check=.TRUE.
269    !
270    !
271      orogfname = 'Relief.nc'
272    !
273      IF ( check ) WRITE(*,*) 'Reading the high resolution orography'
274    !
275      CALL flininfo(orogfname,iml_rel, jml_rel, llm_tmp, ttm_tmp, fid)
276    !
277      ALLOCATE (lat_rel(iml_rel,jml_rel), stat=iret)
278      ALLOCATE (lon_rel(iml_rel,jml_rel), stat=iret)
279      ALLOCATE (relief_hi(iml_rel,jml_rel), stat=iret)
280    !
281      CALL flinopen(orogfname, .FALSE., iml_rel, jml_rel,
282     .llm_tmp, lon_rel, lat_rel, lev, ttm_tmp,
283     .      itau, date, dt, fid)
284    !
285      CALL flinget(fid, 'RELIEF', iml_rel, jml_rel, llm_tmp,
286     . ttm_tmp, 1, 1, relief_hi)
287    !
288      CALL flinclo(fid)
289    !
290    !   In case we have a file which is in degrees we do the transformation
291    !
292      ALLOCATE(lon_rad(iml_rel))
293      ALLOCATE(lon_ini(iml_rel))
294
295      IF ( MAXVAL(lon_rel(:,:)) .GT. 2.0 * ASIN(1.0) ) THEN
296          lon_ini(:) = lon_rel(:,1) * 2.0 * ASIN(1.0) / 180.0
297      ELSE
298          lon_ini(:) = lon_rel(:,1)
299      ENDIF
300
301      ALLOCATE(lat_rad(jml_rel))
302      ALLOCATE(lat_ini(jml_rel))
303
304      IF ( MAXVAL(lat_rel(:,:)) .GT. 2.0 * ASIN(1.0) ) THEN
305          lat_ini(:) = lat_rel(1,:) * 2.0 * ASIN(1.0) / 180.0
306      ELSE
307          lat_ini(:) = lat_rel(1,:)
308      ENDIF
309    !
310    !
311
312      title='RELIEF'
313
314      interbar2 = .FALSE.
315      CALL conf_dat2d(title,iml_rel, jml_rel, lon_ini, lat_ini,
316     . lon_rad, lat_rad, relief_hi , interbar2  )
317
318      IF ( check ) WRITE(*,*) 'Computes all the parameters needed',
319     .' for the gravity wave drag code'
320    !
321    !    Allocate the data we need to put in the interpolated fields
322    !
323    !            RELIEF:  orographie moyenne
324      ALLOCATE(relief(iml,jml))
325    !            zphi :  orographie moyenne
326      ALLOCATE(phis(iml,jml))
327    !             zstd:  deviation standard de l'orographie sous-maille
328      ALLOCATE(zstd(iml,jml))
329    !             zsig:  pente de l'orographie sous-maille
330      ALLOCATE(zsig(iml,jml))
331    !             zgam:  anisotropy de l'orographie sous maille
332      ALLOCATE(zgam(iml,jml))
333    !             zthe:  orientation de l'axe oriente dans la direction
334    !                    de plus grande pente de l'orographie sous maille
335      ALLOCATE(zthe(iml,jml))
336    !             zpic:  hauteur pics de la SSO
337      ALLOCATE(zpic(iml,jml))
338    !             zval:  hauteur vallees de la SSO
339      ALLOCATE(zval(iml,jml))
340    !             masque : Masque terre ocean
341      ALLOCATE(tmp_int(iml,jml))
342      ALLOCATE(masque(iml,jml))
343
344      masque = -99999.
345      if (present(masque_lu)) then
346        masque = masque_lu
347      endif
348    !
349      CALL grid_noro(iml_rel, jml_rel, lon_rad, lat_rad, relief_hi,
350     . iml-1, jml, lon_in, lat_in,
351     . phis, relief, zstd, zsig, zgam, zthe, zpic, zval, masque)
352      phis = phis * 9.81
353    !
354!      masque(:,:) = REAL(tmp_int(:,:))
355    !
356    !  Compute surface roughness
357    !
358      IF ( check ) WRITE(*,*)
359     .'Compute surface roughness induced by the orography'
360    !
361      ALLOCATE(rugo(iml,jml))
362      ALLOCATE(tmp_var(iml-1,jml))
363    !
364      CALL rugsoro(iml_rel, jml_rel, lon_rad, lat_rad, relief_hi,
365     . iml-1, jml, lon_in, lat_in, tmp_var)
366    !
367      DO j = 1, jml
368        DO i = 1, iml-1
369          rugo(i,j) = tmp_var(i,j)
370        ENDDO
371        rugo(iml,j) = tmp_var(1,j)
372      ENDDO
373c
374cc   ***   rugo  n'est pas utilise pour l'instant  ******
375    !
376    !   Build land-sea mask
377    !
378    !
379      RETURN
380    !
381      END SUBROUTINE start_init_orog
382    !
383    !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
384    !
385      SUBROUTINE startget_phys1d(varname, iml, jml, lon_in,
386     .lat_in, nbindex, champ, val_exp ,jml2, lon_in2, lat_in2,interbar)
387    !
388      CHARACTER*(*), INTENT(in) :: varname
389      INTEGER, INTENT(in) :: iml, jml, nbindex, jml2
390      REAL, INTENT(in) :: lon_in(iml), lat_in(jml)
391      REAL, INTENT(in) :: lon_in2(iml), lat_in2(jml2)
392      REAL, INTENT(inout) :: champ(nbindex)
393      REAL, INTENT(in) :: val_exp
394      LOGICAL interbar
395    !
396    !
397    !   This routine only works if the variable does not exist or is constant
398    !
399      IF ( MINVAL(champ(:)).EQ.MAXVAL(champ(:)) .AND.
400     .MINVAL(champ(:)).EQ.val_exp ) THEN
401          SELECTCASE(varname)
402            CASE ('tsol')
403              IF ( .NOT.ALLOCATED(tsol)) THEN
404                CALL start_init_phys( iml, jml, lon_in, lat_in,
405     .              jml2, lon_in2, lat_in2, interbar )
406              ENDIF
407              IF ( SIZE(tsol) .NE. SIZE(lon_in)*SIZE(lat_in) ) THEN
408                WRITE(*,*)
409     . 'STARTVAR module has been initialized to the wrong size'
410                 STOP
411              ENDIF
412              CALL gr_dyn_fi(1, iml, jml, nbindex, tsol, champ)
413            CASE ('qsol')
414              IF ( .NOT.ALLOCATED(qsol)) THEN
415                CALL start_init_phys( iml, jml, lon_in, lat_in,
416     .              jml2, lon_in2,lat_in2 , interbar )
417              ENDIF
418              IF ( SIZE(qsol) .NE. SIZE(lon_in)*SIZE(lat_in) ) THEN
419                WRITE(*,*)
420     . 'STARTVAR module has been initialized to the wrong size'
421                STOP
422              ENDIF
423              CALL gr_dyn_fi(1, iml, jml, nbindex, qsol, champ)
424            CASE ('psol')
425              IF ( .NOT.ALLOCATED(psol_dyn)) THEN
426                CALL start_init_dyn( iml, jml, lon_in, lat_in,
427     .              jml2, lon_in2,lat_in2 , interbar )
428              ENDIF
429              IF (SIZE(psol_dyn) .NE. SIZE(lon_in)*SIZE(lat_in)) THEN
430                WRITE(*,*)
431     . 'STARTVAR module has been initialized to the wrong size'
432                STOP
433              ENDIF
434              CALL gr_dyn_fi(1, iml, jml, nbindex, psol_dyn, champ)
435            CASE ('zmea')
436              IF ( .NOT.ALLOCATED(relief)) THEN
437                CALL start_init_orog( iml, jml, lon_in, lat_in,
438     .            jml2, lon_in2,lat_in2 , interbar )
439              ENDIF
440              IF ( SIZE(relief) .NE. SIZE(lon_in)*SIZE(lat_in) ) THEN
441                WRITE(*,*)
442     . 'STARTVAR module has been initialized to the wrong size'
443                 STOP
444              ENDIF
445              CALL gr_dyn_fi(1, iml, jml, nbindex, relief, champ)
446            CASE ('zstd')
447              IF ( .NOT.ALLOCATED(zstd)) THEN
448                CALL start_init_orog( iml, jml, lon_in, lat_in,
449     .              jml2, lon_in2,lat_in2 , interbar )
450              ENDIF
451              IF ( SIZE(zstd) .NE. SIZE(lon_in)*SIZE(lat_in) ) THEN
452                WRITE(*,*)
453     . 'STARTVAR module has been initialized to the wrong size'
454                 STOP
455              ENDIF
456              CALL gr_dyn_fi(1, iml, jml, nbindex,zstd, champ)
457            CASE ('zsig')
458              IF ( .NOT.ALLOCATED(zsig)) THEN
459                CALL start_init_orog( iml, jml, lon_in, lat_in,
460     .               jml2, lon_in2,lat_in2 , interbar )
461              ENDIF
462              IF ( SIZE(zsig) .NE. SIZE(lon_in)*SIZE(lat_in) ) THEN
463                WRITE(*,*)
464     . 'STARTVAR module has been initialized to the wrong size'
465                 STOP
466              ENDIF
467              CALL gr_dyn_fi(1, iml, jml, nbindex,zsig, champ)
468            CASE ('zgam')
469              IF ( .NOT.ALLOCATED(zgam)) THEN
470                CALL start_init_orog( iml, jml, lon_in, lat_in,
471     .            jml2, lon_in2,lat_in2 , interbar )
472              ENDIF
473              IF ( SIZE(zgam) .NE. SIZE(lon_in)*SIZE(lat_in) ) THEN
474                WRITE(*,*)
475     . 'STARTVAR module has been initialized to the wrong size'
476                 STOP
477              ENDIF
478              CALL gr_dyn_fi(1, iml, jml, nbindex,zgam, champ)
479            CASE ('zthe')
480              IF ( .NOT.ALLOCATED(zthe)) THEN
481                CALL start_init_orog( iml, jml, lon_in, lat_in,
482     .            jml2, lon_in2,lat_in2 , interbar )
483              ENDIF
484              IF ( SIZE(zthe) .NE. SIZE(lon_in)*SIZE(lat_in) ) THEN
485                WRITE(*,*)
486     . 'STARTVAR module has been initialized to the wrong size'
487                 STOP
488              ENDIF
489              CALL gr_dyn_fi(1, iml, jml, nbindex,zthe, champ)
490            CASE ('zpic')
491              IF ( .NOT.ALLOCATED(zpic)) THEN
492                CALL start_init_orog( iml, jml, lon_in, lat_in,
493     .            jml2, lon_in2,lat_in2 , interbar )
494              ENDIF
495              IF ( SIZE(zpic) .NE. SIZE(lon_in)*SIZE(lat_in) ) THEN
496                WRITE(*,*)
497     . 'STARTVAR module has been initialized to the wrong size'
498                 STOP
499              ENDIF
500              CALL gr_dyn_fi(1, iml, jml, nbindex,zpic, champ)
501            CASE ('zval')
502              IF ( .NOT.ALLOCATED(zval)) THEN
503                CALL start_init_orog( iml, jml, lon_in, lat_in,
504     .            jml2, lon_in2,lat_in2 , interbar )
505              ENDIF
506              IF ( SIZE(zval) .NE. SIZE(lon_in)*SIZE(lat_in) ) THEN
507                WRITE(*,*)
508     . 'STARTVAR module has been initialized to the wrong size'
509                 STOP
510              ENDIF
511              CALL gr_dyn_fi(1, iml, jml, nbindex,zval, champ)
512            CASE ('rads')
513                  champ(:) = 0.0
514            CASE ('snow')
515                  champ(:) = 0.0
516cIM "slab" ocean
517            CASE ('tslab')
518                   champ(:) = 0.0
519            CASE ('seaice')
520                  champ(:) = 0.0
521            CASE ('rugmer')
522                  champ(:) = 0.001
523            CASE ('agsno')
524                  champ(:) = 50.0
525            CASE DEFAULT
526              WRITE(*,*) 'startget_phys1d'
527              WRITE(*,*) 'No rule is present to extract variable  ',
528     . varname(:LEN_TRIM(varname)),' from any data set'
529              STOP
530          END SELECT
531      ELSE
532        !
533        ! If we see tsol we catch it as we may need it for a 3D interpolation
534        !
535        SELECTCASE(varname)
536          CASE ('tsol')
537            IF ( .NOT.ALLOCATED(tsol)) THEN
538              ALLOCATE(tsol(SIZE(lon_in),SIZE(lat_in) ))
539            ENDIF
540            CALL gr_fi_dyn(1, iml, jml, nbindex, champ, tsol)
541        END SELECT
542      ENDIF
543      END SUBROUTINE startget_phys1d
544    !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
545    !
546      SUBROUTINE start_init_phys( iml, jml, lon_in, lat_in, jml2,
547     .                 lon_in2, lat_in2 , interbar )
548
549      use inter_barxy_m, only: inter_barxy
550    !
551      INTEGER, INTENT(in) :: iml, jml ,jml2
552      REAL, INTENT(in) :: lon_in(iml), lat_in(jml)
553      REAL, INTENT(in) :: lon_in2(iml), lat_in2(jml2)
554      LOGICAL interbar
555    !
556    !  LOCAL
557    !
558!ac     REAL :: lev(1), date, dt
559      REAL :: date, dt
560      REAL, DIMENSION(:), ALLOCATABLE :: levphys_ini
561!ac
562      INTEGER :: itau(1)
563      INTEGER ::  llm_tmp, ttm_tmp
564      INTEGER :: i, j
565    !
566      CHARACTER*25 title
567      CHARACTER*120 :: physfname
568      LOGICAL :: check=.TRUE.
569    !
570      REAL, ALLOCATABLE :: lon_rad(:), lat_rad(:)
571      REAL, ALLOCATABLE :: lon_ini(:), lat_ini(:)
572      REAL, ALLOCATABLE :: var_ana(:,:), tmp_var(:,:)
573    !
574      physfname = 'ECPHY.nc'
575    !
576      IF ( check ) WRITE(*,*) 'Opening the surface analysis'
577    !
578      CALL flininfo(physfname, iml_phys, jml_phys, llm_tmp,
579     . ttm_tmp, fid_phys)
580    !
581      ALLOCATE (lat_phys(iml_phys,jml_phys))
582      ALLOCATE (lon_phys(iml_phys,jml_phys))
583!ac
584      ALLOCATE (levphys_ini(llm_tmp))
585    !
586!      CALL flinopen(physfname, .FALSE., iml_phys, jml_phys,
587!     . llm_tmp, lon_phys, lat_phys, lev, ttm_tmp,
588!     . itau, date, dt, fid_phys)
589    !
590      CALL flinopen(physfname, .FALSE., iml_phys, jml_phys,
591     . llm_tmp, lon_phys, lat_phys, levphys_ini, ttm_tmp,
592     . itau, date, dt, fid_phys)
593    !
594      DEALLOCATE (levphys_ini)
595!ac
596    !
597    ! Allocate the space we will need to get the data out of this file
598    !
599      ALLOCATE(var_ana(iml_phys, jml_phys))
600    !
601    !   In case we have a file which is in degrees we do the transformation
602    !
603      ALLOCATE(lon_rad(iml_phys))
604      ALLOCATE(lon_ini(iml_phys))
605
606      IF ( MAXVAL(lon_phys(:,:)) .GT. 2.0 * ASIN(1.0) ) THEN
607          lon_ini(:) = lon_phys(:,1) * 2.0 * ASIN(1.0) / 180.0
608      ELSE
609          lon_ini(:) = lon_phys(:,1)
610      ENDIF
611
612      ALLOCATE(lat_rad(jml_phys))
613      ALLOCATE(lat_ini(jml_phys))
614
615      IF ( MAXVAL(lat_phys(:,:)) .GT. 2.0 * ASIN(1.0) ) THEN
616          lat_ini(:) = lat_phys(1,:) * 2.0 * ASIN(1.0) / 180.0
617      ELSE
618          lat_ini(:) = lat_phys(1,:)
619      ENDIF
620
621
622    !
623    !   We get the two standard varibales
624    !   Surface temperature
625    !
626      ALLOCATE(tsol(iml,jml))
627      ALLOCATE(tmp_var(iml-1,jml))
628    !
629    !
630
631      CALL flinget(fid_phys, 'ST', iml_phys, jml_phys,
632     .llm_tmp, ttm_tmp, 1, 1, var_ana)
633
634      title='ST'
635      CALL conf_dat2d(title,iml_phys, jml_phys, lon_ini, lat_ini,
636     . lon_rad, lat_rad, var_ana , interbar  )
637
638      IF ( interbar )   THEN
639        WRITE(6,*) '-------------------------------------------------',
640     ,'--------------'
641        WRITE(6,*) '$$$ Utilisation de l interpolation barycentrique ',
642     , ' pour  ST $$$ '
643        WRITE(6,*) '-------------------------------------------------',
644     ,'--------------'
645        CALL inter_barxy(lon_rad, lat_rad(:jml_phys -1), var_ana,
646     $       lon_in2(:iml-1), lat_in2(:jml-1), tmp_var)
647      ELSE
648        CALL grille_m(iml_phys, jml_phys, lon_rad, lat_rad,
649     .    var_ana, iml-1, jml, lon_in, lat_in, tmp_var     )
650      ENDIF
651
652      CALL gr_int_dyn(tmp_var, tsol, iml-1, jml)
653    !
654    ! Soil moisture
655    !
656      ALLOCATE(qsol(iml,jml))
657      CALL flinget(fid_phys, 'CDSW', iml_phys, jml_phys,
658     . llm_tmp, ttm_tmp, 1, 1, var_ana)
659
660      title='CDSW'
661      CALL conf_dat2d(title,iml_phys, jml_phys, lon_ini, lat_ini,
662     . lon_rad, lat_rad, var_ana, interbar  )
663
664      IF ( interbar )   THEN
665        WRITE(6,*) '-------------------------------------------------',
666     ,'--------------'
667        WRITE(6,*) '$$$ Utilisation de l interpolation barycentrique ',
668     , ' pour  CDSW $$$ '
669        WRITE(6,*) '-------------------------------------------------',
670     ,'--------------'
671        CALL inter_barxy(lon_rad, lat_rad(:jml_phys -1), var_ana,
672     $       lon_in2(:iml-1), lat_in2(:jml-1), tmp_var)
673      ELSE
674        CALL grille_m(iml_phys, jml_phys, lon_rad, lat_rad,
675     .    var_ana, iml-1, jml, lon_in, lat_in, tmp_var     )
676      ENDIF
677c
678        CALL gr_int_dyn(tmp_var, qsol, iml-1, jml)
679    !
680       CALL flinclo(fid_phys)
681    !
682      END SUBROUTINE start_init_phys
683    !
684    !
685    !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
686    !
687    !
688      SUBROUTINE startget_dyn(varname, lon_in, lat_in,
689     . pls, workvar, champ, val_exp, lon_in2, lat_in2 ,
690     ,  interbar )
691
692      use assert_eq_m, only: assert_eq
693    !
694    !   ARGUMENTS
695    !
696      CHARACTER(len=*), INTENT(in) :: varname
697      REAL, INTENT(in) :: lon_in(:) ! dim(iml)
698      REAL, INTENT(in) :: lat_in(:) ! dim(jml)
699      REAL, INTENT(in) :: lon_in2(:) ! dim(iml)
700      REAL, INTENT(in) :: lat_in2(:) ! dim(jml2)
701      REAL, INTENT(in) :: pls(:, :, :) ! dim(iml, jml, lml)
702      REAL, INTENT(in) :: workvar(:, :, :) ! dim(iml, jml, lml)
703      REAL, INTENT(inout) :: champ(:, :, :) ! dim(iml, jml, lml)
704      REAL, INTENT(in) :: val_exp
705      LOGICAL interbar
706    !
707    !    LOCAL
708    !
709      INTEGER :: il, ij, ii, iml, jml, lml, jml2
710      REAL :: xppn, xpps
711    !
712#include "dimensions.h"
713#include "paramet.h"
714#include "comgeom2.h"
715#include "comconst.h"
716    !
717    !   This routine only works if the variable does not exist or is constant
718    !
719C     -----------------------------
720
721      iml = assert_eq((/size(lon_in), size(pls, 1), size(workvar, 1),
722     $     size(champ, 1), size(lon_in2)/), "startget_dyn iml")
723      jml = assert_eq(size(lat_in), size(pls, 2), size(workvar, 2),
724     $     size(champ, 2), "startget_dyn jml")
725      lml = assert_eq(size(pls, 3), size(workvar, 3), size(champ, 3),
726     $     "startget_dyn lml")
727      jml2 = size(lat_in2)
728
729      IF ( MINVAL(champ(:,:,:)).EQ.MAXVAL(champ(:,:,:)) .AND.
730     . MINVAL(champ(:,:,:)).EQ.val_exp ) THEN
731        !
732        SELECTCASE(varname)
733          CASE ('u')
734            IF ( .NOT.ALLOCATED(psol_dyn)) THEN
735              CALL start_init_dyn( iml, jml, lon_in, lat_in, jml2 ,
736     .          lon_in2,lat_in2 , interbar )
737            ENDIF
738            CALL start_inter_3d('U', iml, jml, lml, lon_in,
739     .       lat_in, jml2, lon_in2, lat_in2,  pls, champ,interbar )
740            DO il=1,lml
741              DO ij=1,jml
742                DO ii=1,iml-1
743                  champ(ii,ij,il) = champ(ii,ij,il) * cu(ii,ij)
744                ENDDO
745                champ(iml,ij, il) = champ(1,ij, il)
746              ENDDO
747            ENDDO
748          CASE ('v')
749            IF ( .NOT.ALLOCATED(psol_dyn)) THEN
750              CALL start_init_dyn( iml, jml, lon_in, lat_in , jml2,
751     .           lon_in2, lat_in2 , interbar )
752            ENDIF
753            CALL start_inter_3d('V', iml, jml, lml, lon_in,
754     .       lat_in, jml2, lon_in2, lat_in2,  pls, champ, interbar )
755            DO il=1,lml
756              DO ij=1,jml
757                DO ii=1,iml-1
758                  champ(ii,ij,il) = champ(ii,ij,il) * cv(ii,ij)
759                ENDDO
760                champ(iml,ij, il) = champ(1,ij, il)
761              ENDDO
762            ENDDO
763          CASE ('t')
764            IF ( .NOT.ALLOCATED(psol_dyn)) THEN
765              CALL start_init_dyn( iml, jml, lon_in, lat_in, jml2 ,
766     .           lon_in2, lat_in2 ,interbar )
767            ENDIF
768            CALL start_inter_3d('TEMP', iml, jml, lml, lon_in,
769     .       lat_in, jml2, lon_in2, lat_in2,  pls, champ, interbar )
770 
771          CASE ('tpot')
772            IF ( .NOT.ALLOCATED(psol_dyn)) THEN
773              CALL start_init_dyn( iml, jml, lon_in, lat_in , jml2 ,
774     .            lon_in2, lat_in2 , interbar )
775            ENDIF
776            CALL start_inter_3d('TEMP', iml, jml, lml, lon_in,
777     .       lat_in, jml2, lon_in2, lat_in2,  pls, champ, interbar )
778            IF ( MINVAL(workvar(:,:,:)) .NE. MAXVAL(workvar(:,:,:)) )
779     .                                    THEN
780              DO il=1,lml
781                DO ij=1,jml
782                  DO ii=1,iml-1
783                    champ(ii,ij,il) = champ(ii,ij,il) * cpp
784     .                                 / workvar(ii,ij,il)
785                  ENDDO
786                  champ(iml,ij,il) = champ(1,ij,il)
787                ENDDO
788              ENDDO
789              DO il=1,lml
790                xppn = SUM(aire(:,1)*champ(:,1,il))/apoln
791                xpps = SUM(aire(:,jml)*champ(:,jml,il))/apols
792                champ(:,1,il) = xppn
793                champ(:,jml,il) = xpps
794              ENDDO
795            ELSE
796              WRITE(*,*)'Could not compute potential temperature as the'
797              WRITE(*,*)'Exner function is missing or constant.'
798              STOP
799            ENDIF
800          CASE ('q')
801            IF ( .NOT.ALLOCATED(psol_dyn)) THEN
802              CALL start_init_dyn( iml, jml, lon_in, lat_in, jml2 ,
803     .           lon_in2, lat_in2 , interbar )
804            ENDIF
805            CALL start_inter_3d('R', iml, jml, lml, lon_in, lat_in,
806     .        jml2, lon_in2, lat_in2,  pls, champ, interbar )
807            IF ( MINVAL(workvar(:,:,:)) .NE. MAXVAL(workvar(:,:,:)) )
808     .                                     THEN
809              DO il=1,lml
810                DO ij=1,jml
811                  DO ii=1,iml-1
812                    champ(ii,ij,il) = 0.01 * champ(ii,ij,il) *
813     .                                       workvar(ii,ij,il)
814                  ENDDO
815                  champ(iml,ij,il) = champ(1,ij,il)
816                ENDDO
817              ENDDO
818              WHERE ( champ .LT. 0.) champ = 1.0E-10
819              DO il=1,lml
820                xppn = SUM(aire(:,1)*champ(:,1,il))/apoln
821                xpps = SUM(aire(:,jml)*champ(:,jml,il))/apols
822                champ(:,1,il) = xppn
823                champ(:,jml,il) = xpps
824              ENDDO
825            ELSE
826              WRITE(*,*)'Could not compute specific humidity as the'
827              WRITE(*,*)'saturated humidity is missing or constant.'
828              STOP
829            ENDIF
830          CASE DEFAULT
831            WRITE(*,*) 'startget_dyn'
832            WRITE(*,*) 'No rule is present to extract variable  ',
833     . varname(:LEN_TRIM(varname)),' from any data set'
834            STOP
835          END SELECT
836      ENDIF
837      END SUBROUTINE startget_dyn
838    !
839    !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
840    !
841      SUBROUTINE start_init_dyn( iml, jml, lon_in, lat_in,jml2,lon_in2 ,
842     ,             lat_in2 , interbar )
843    !
844      use inter_barxy_m, only: inter_barxy
845
846      INTEGER, INTENT(in) :: iml, jml, jml2
847      REAL, INTENT(in) :: lon_in(iml), lat_in(jml)
848      REAL, INTENT(in) :: lon_in2(iml), lat_in2(jml2)
849      LOGICAL interbar
850    !
851    !  LOCAL
852    !
853      REAL :: lev(1), date, dt
854      INTEGER :: itau(1)
855      INTEGER :: i, j
856      integer :: iret
857    !
858      CHARACTER*120 :: physfname
859      LOGICAL :: check=.TRUE.
860    !
861      REAL, ALLOCATABLE :: lon_rad(:), lat_rad(:)
862      REAL, ALLOCATABLE :: lon_ini(:), lat_ini(:)
863      REAL, ALLOCATABLE :: var_ana(:,:), tmp_var(:,:), z(:,:)
864      REAL, ALLOCATABLE :: xppn(:), xpps(:)
865      LOGICAL :: allo
866    !
867    !
868#include "dimensions.h"
869#include "paramet.h"
870#include "comgeom2.h"
871
872      CHARACTER*25 title
873
874    !
875      physfname = 'ECDYN.nc'
876    !
877      IF ( check ) WRITE(*,*) 'Opening the surface analysis'
878    !
879      CALL flininfo(physfname, iml_dyn, jml_dyn, llm_dyn,
880     .                            ttm_dyn, fid_dyn)
881      IF ( check ) WRITE(*,*) 'Values read: ', iml_dyn, jml_dyn,
882     .                                         llm_dyn, ttm_dyn
883    !
884      ALLOCATE (lat_dyn(iml_dyn,jml_dyn), stat=iret)
885      ALLOCATE (lon_dyn(iml_dyn,jml_dyn), stat=iret)
886      ALLOCATE (levdyn_ini(llm_dyn), stat=iret)
887    !
888      CALL flinopen(physfname, .FALSE., iml_dyn, jml_dyn, llm_dyn,
889     . lon_dyn, lat_dyn, levdyn_ini, ttm_dyn,
890     . itau, date, dt, fid_dyn)
891    !
892
893      allo = allocated (var_ana)
894      if (allo) then
895        DEALLOCATE(var_ana, stat=iret)
896      endif
897      ALLOCATE(var_ana(iml_dyn, jml_dyn), stat=iret)
898
899      allo = allocated (lon_rad)
900      if (allo) then
901        DEALLOCATE(lon_rad, stat=iret)
902      endif
903
904      ALLOCATE(lon_rad(iml_dyn), stat=iret)
905      ALLOCATE(lon_ini(iml_dyn))
906       
907      IF ( MAXVAL(lon_dyn(:,:)) .GT. 2.0 * ASIN(1.0) ) THEN
908          lon_ini(:) = lon_dyn(:,1) * 2.0 * ASIN(1.0) / 180.0
909      ELSE
910          lon_ini(:) = lon_dyn(:,1)
911      ENDIF
912
913      ALLOCATE(lat_rad(jml_dyn))
914      ALLOCATE(lat_ini(jml_dyn))
915
916      IF ( MAXVAL(lat_dyn(:,:)) .GT. 2.0 * ASIN(1.0) ) THEN
917          lat_ini(:) = lat_dyn(1,:) * 2.0 * ASIN(1.0) / 180.0
918      ELSE
919          lat_ini(:) = lat_dyn(1,:)
920      ENDIF
921    !
922
923
924      ALLOCATE(z(iml, jml))
925      ALLOCATE(tmp_var(iml-1,jml))
926    !
927      CALL flinget(fid_dyn, 'Z', iml_dyn, jml_dyn, 0, ttm_dyn,
928     .              1, 1, var_ana)
929c
930      title='Z'
931      CALL conf_dat2d( title,iml_dyn, jml_dyn,lon_ini, lat_ini,
932     . lon_rad, lat_rad, var_ana, interbar  )
933c
934      IF ( interbar )   THEN
935        WRITE(6,*) '-------------------------------------------------',
936     ,'--------------'
937        WRITE(6,*) '$$$ Utilisation de l interpolation barycentrique ',
938     , ' pour  Z  $$$ '
939        WRITE(6,*) '-------------------------------------------------',
940     ,'--------------'
941        CALL inter_barxy(lon_rad, lat_rad(:jml_dyn -1), var_ana,
942     $       lon_in2(:iml-1), lat_in2(:jml-1), tmp_var)
943      ELSE
944        CALL grille_m(iml_dyn, jml_dyn , lon_rad, lat_rad, var_ana,
945     .               iml-1, jml, lon_in, lat_in, tmp_var)
946      ENDIF
947
948      CALL gr_int_dyn(tmp_var, z, iml-1, jml)
949    !
950      ALLOCATE(psol_dyn(iml, jml))
951    !
952      CALL flinget(fid_dyn, 'SP', iml_dyn, jml_dyn, 0, ttm_dyn,
953     .              1, 1, var_ana)
954
955       title='SP'
956      CALL conf_dat2d( title,iml_dyn, jml_dyn,lon_ini, lat_ini,
957     . lon_rad, lat_rad, var_ana, interbar  )
958
959      IF ( interbar )   THEN
960        WRITE(6,*) '-------------------------------------------------',
961     ,'--------------'
962        WRITE(6,*) '$$$ Utilisation de l interpolation barycentrique ',
963     , ' pour  SP  $$$ '
964        WRITE(6,*) '-------------------------------------------------',
965     ,'--------------'
966        CALL inter_barxy(lon_rad, lat_rad(:jml_dyn -1), var_ana,
967     $       lon_in2(:iml-1), lat_in2(:jml-1), tmp_var)
968      ELSE
969        CALL grille_m(iml_dyn, jml_dyn , lon_rad, lat_rad, var_ana,
970     .             iml-1, jml, lon_in, lat_in, tmp_var  )
971      ENDIF
972
973      CALL gr_int_dyn(tmp_var, psol_dyn, iml-1, jml)
974    !
975      IF ( .NOT.ALLOCATED(tsol)) THEN
976    !   These variables may have been allocated by the need to
977    !   create a start field for them or by the varibale
978    !   coming out of the restart file. In case we dor have it we will initialize it.
979    !
980        CALL start_init_phys( iml, jml, lon_in, lat_in,jml2,lon_in2,
981     .                 lat_in2 , interbar )
982      ELSE
983        IF ( SIZE(tsol) .NE. SIZE(psol_dyn) ) THEN
984        WRITE(*,*) 'start_init_dyn :'
985        WRITE(*,*) 'The temperature field we have does not ',
986     .             'have the right size'
987        STOP
988      ENDIF
989      ENDIF
990      IF ( .NOT.ALLOCATED(phis)) THEN
991            !
992            !    These variables may have been allocated by the need to create a start field for them or by the varibale
993            !     coming out of the restart file. In case we dor have it we will initialize it.
994            !
995        CALL start_init_orog( iml, jml, lon_in, lat_in, jml2, lon_in2 ,
996     .      lat_in2 , interbar )
997            !
998      ELSE
999            !
1000          IF (SIZE(phis) .NE. SIZE(psol_dyn)) THEN
1001                !
1002              WRITE(*,*) 'start_init_dyn :'
1003              WRITE(*,*) 'The orography field we have does not ',
1004     .                   ' have the right size'
1005              STOP
1006          ENDIF
1007            !
1008      ENDIF
1009    !
1010    !     PSOL is computed in Pascals
1011    !
1012    !
1013      DO j = 1, jml
1014        DO i = 1, iml-1
1015          psol_dyn(i,j) = psol_dyn(i,j)*(1.0+(z(i,j)-phis(i,j))
1016     .                    /287.0/tsol(i,j))
1017        ENDDO
1018        psol_dyn(iml,j) = psol_dyn(1,j)
1019      ENDDO
1020    !
1021    !
1022      ALLOCATE(xppn(iml-1))
1023      ALLOCATE(xpps(iml-1))
1024    !
1025      DO  i   = 1, iml-1
1026        xppn(i) = aire( i,1) * psol_dyn( i,1)
1027        xpps(i) = aire( i,jml) * psol_dyn( i,jml)
1028      ENDDO
1029    !
1030      DO i   = 1, iml
1031        psol_dyn(i,1    )  = SUM(xppn)/apoln
1032        psol_dyn(i,jml)  = SUM(xpps)/apols
1033      ENDDO
1034    !
1035      RETURN
1036    !
1037      END SUBROUTINE start_init_dyn
1038    !
1039    !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1040    !
1041      SUBROUTINE start_inter_3d(varname, iml, jml, lml, lon_in,
1042     .      lat_in, jml2, lon_in2, lat_in2, pls_in, var3d, interbar )
1043    !
1044    !    This subroutine gets a variables from a 3D file and does the interpolations needed
1045    !
1046    !
1047      use inter_barxy_m, only: inter_barxy
1048
1049    !    ARGUMENTS
1050    !
1051      CHARACTER*(*) :: varname
1052      INTEGER :: iml, jml, lml, jml2
1053      REAL :: lon_in(iml), lat_in(jml), pls_in(iml, jml, lml)
1054      REAL :: lon_in2(iml) , lat_in2(jml2)
1055      REAL :: var3d(iml, jml, lml)
1056      LOGICAL interbar
1057      real chmin,chmax
1058    !
1059    !  LOCAL
1060    !
1061      CHARACTER*25 title
1062      INTEGER :: ii, ij, il, jsort,i,j,l
1063      REAL :: bx, by
1064      REAL, ALLOCATABLE :: lon_rad(:), lat_rad(:)
1065      REAL, ALLOCATABLE :: lon_ini(:), lat_ini(:) , lev_dyn(:)
1066      REAL, ALLOCATABLE :: var_tmp2d(:,:), var_tmp3d(:,:,:)
1067      REAL, ALLOCATABLE :: ax(:), ay(:), yder(:)
1068!       REAL, ALLOCATABLE :: varrr(:,:,:)
1069      INTEGER, ALLOCATABLE :: lind(:)
1070    !
1071      LOGICAL :: check = .TRUE.
1072    !
1073      IF ( .NOT. ALLOCATED(var_ana3d)) THEN
1074          ALLOCATE(var_ana3d(iml_dyn, jml_dyn, llm_dyn))
1075      ENDIF
1076!          ALLOCATE(varrr(iml_dyn, jml_dyn, llm_dyn))
1077    !
1078    !
1079      IF ( check) WRITE(*,*) 'Going into flinget to extract the 3D ',
1080     .  ' field.', fid_dyn
1081      IF ( check) WRITE(*,*) fid_dyn, varname, iml_dyn, jml_dyn,
1082     .                        llm_dyn,ttm_dyn
1083    !
1084      CALL flinget(fid_dyn, varname, iml_dyn, jml_dyn, llm_dyn,
1085     . ttm_dyn, 1, 1, var_ana3d)
1086    !
1087      IF ( check) WRITE(*,*) 'Allocating space for the interpolation',
1088     . iml, jml, llm_dyn
1089    !
1090      ALLOCATE(lon_rad(iml_dyn))
1091      ALLOCATE(lon_ini(iml_dyn))
1092
1093      IF ( MAXVAL(lon_dyn(:,:)) .GT. 2.0 * ASIN(1.0) ) THEN
1094          lon_ini(:) = lon_dyn(:,1) * 2.0 * ASIN(1.0) / 180.0
1095      ELSE
1096          lon_ini(:) = lon_dyn(:,1)
1097      ENDIF
1098
1099      ALLOCATE(lat_rad(jml_dyn))
1100      ALLOCATE(lat_ini(jml_dyn))
1101
1102      ALLOCATE(lev_dyn(llm_dyn))
1103
1104      IF ( MAXVAL(lat_dyn(:,:)) .GT. 2.0 * ASIN(1.0) ) THEN
1105          lat_ini(:) = lat_dyn(1,:) * 2.0 * ASIN(1.0) / 180.0
1106      ELSE
1107          lat_ini(:) = lat_dyn(1,:)
1108      ENDIF
1109    !
1110
1111      CALL conf_dat3d ( varname,iml_dyn, jml_dyn, llm_dyn, lon_ini,
1112     . lat_ini, levdyn_ini, lon_rad, lat_rad, lev_dyn, var_ana3d  ,
1113     ,  interbar                                                   )
1114
1115      ALLOCATE(var_tmp2d(iml-1, jml))
1116      ALLOCATE(var_tmp3d(iml, jml, llm_dyn))
1117      ALLOCATE(ax(llm_dyn))
1118      ALLOCATE(ay(llm_dyn))
1119      ALLOCATE(yder(llm_dyn))
1120      ALLOCATE(lind(llm_dyn))
1121    !
1122 
1123      DO il=1,llm_dyn
1124        !
1125      IF( interbar )  THEN
1126       IF( il.EQ.1 )  THEN
1127        WRITE(6,*) '-------------------------------------------------',
1128     ,'--------------'
1129        WRITE(6,*) '$$$ Utilisation de l interpolation barycentrique ',
1130     , ' pour ', varname
1131        WRITE(6,*) '-------------------------------------------------',
1132     ,'--------------'
1133       ENDIF
1134       CALL inter_barxy(lon_rad, lat_rad(:jml_dyn -1),
1135     $      var_ana3d(:,:,il), lon_in2(:iml-1), lat_in2, var_tmp2d)
1136      ELSE
1137       CALL grille_m(iml_dyn, jml_dyn, lon_rad, lat_rad,
1138     .  var_ana3d(:,:,il), iml-1, jml, lon_in, lat_in, var_tmp2d )
1139      ENDIF
1140        !
1141        CALL gr_int_dyn(var_tmp2d, var_tmp3d(:,:,il), iml-1, jml)
1142        !
1143       ENDDO
1144       !
1145          DO il=1,llm_dyn
1146            lind(il) = llm_dyn-il+1
1147          ENDDO
1148    !
1149c
1150c  ... Pour l'interpolation verticale ,on interpole du haut de l'atmosphere
1151c                    vers  le  sol  ...
1152c
1153      DO ij=1,jml
1154        DO ii=1,iml-1
1155          !
1156          ax(:) = lev_dyn(lind(:))
1157          ay(:) = var_tmp3d(ii, ij, lind(:))
1158          !
1159         
1160          CALL SPLINE(ax, ay, llm_dyn, 1.e30, 1.e30, yder)
1161          !
1162          DO il=1,lml
1163            bx = pls_in(ii, ij, il)
1164            CALL SPLINT(ax, ay, yder, llm_dyn, bx, by)
1165            var3d(ii, ij, il) = by
1166          ENDDO
1167          !
1168        ENDDO
1169        var3d(iml, ij, :) = var3d(1, ij, :)
1170      ENDDO
1171
1172      do il=1,lml
1173        call minmax(iml*jml,var3d(1,1,il),chmin,chmax)
1174      SELECTCASE(varname)
1175       CASE('U')
1176          WRITE(*,*) ' U  min max l ',il,chmin,chmax
1177       CASE('V')
1178          WRITE(*,*) ' V  min max l ',il,chmin,chmax
1179       CASE('TEMP')
1180          WRITE(*,*) ' TEMP  min max l ',il,chmin,chmax
1181       CASE('R')
1182          WRITE(*,*) ' R  min max l ',il,chmin,chmax
1183      END SELECT
1184           enddo
1185
1186      DEALLOCATE(lon_rad)
1187      DEALLOCATE(lon_ini)
1188      DEALLOCATE(lat_rad)
1189      DEALLOCATE(lat_ini)
1190      DEALLOCATE(lev_dyn)
1191      DEALLOCATE(var_tmp2d)
1192      DEALLOCATE(var_tmp3d)
1193      DEALLOCATE(ax)
1194      DEALLOCATE(ay)
1195      DEALLOCATE(yder)
1196      DEALLOCATE(lind)
1197
1198    !
1199      RETURN
1200    !
1201      END SUBROUTINE start_inter_3d
1202    !
1203#endif
1204! of #ifdef CPP_EARTH
1205      END MODULE startvar
Note: See TracBrowser for help on using the repository browser.