source: LMDZ.3.3/branches/rel-LF/libf/dyn3d/startvar.F @ 5444

Last change on this file since 5444 was 488, checked in by lmdzadmin, 21 years ago

Probleme de dimensionnement sur lev (A. Caubel)
LF

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