Philippe 23/02/2016: lfi2cdf: modif for JPHEXT/=1
[MNH-git_open_source-lfs.git] / LIBTOOLS / tools / diachro / src / EXTRACTDIA / concat_time_diafile.f90
1       PROGRAM  EXTRACTDIA
2 !     ###################
3 !
4 !!****  *EXTRACTDIA* -  lecture d'enregistrements dans fichier diachronique,
5 !                         traitement,
6 !                         ecriture (11 types de format de fichier possibles)
7 !! 
8 !!
9 !!    PURPOSE
10 !!    -------
11
12 !
13 !!**  METHOD
14 !!    ------
15 !!      
16 !     Lecture en entree:
17 !       d'une liste de fichiers diachroniques 
18 !       du format de sortie
19 !       d'une liste de champs a traiter pour chaque fichier diachronique
20 !       d'un zoom selon toutes les directions inclu dans le champ a traiter
21 !         ( seul le zoom selon i,j,k est possible pour le format DIAC)
22 !
23 !     Ecriture en sortie:
24 !       d'un fichier  au format fonction de TYPEOUT c.a.d 
25 !         DIAC= type diachro (un seul fichier contenant toutes
26 !                                       les variables selectionnées)
27 !         LLHV= lon lat alt val (un seul fichier contenant toutes
28 !                                       les variables selectionnées) 
29 !         llhv= lat lon alt val (un seul fichier contenant toutes
30 !                                       les variables selectionnées) 
31 !         ll ou LL zv lon lat  niveau Z val
32 !
33 !         ll ou LL pv lon lat  niveau P val
34 !         FREE= format libre a choisir par l utilisateur (un fichier par variable)
35 !         KCDL ou ZCDL ou PCDL= format CDL (à convertir en netcdf via "tonetcdf")
36 !                               (un seul fichier contenant toutes
37 !                                       les variables selectionnées)
38 !           KCDL si les niveaux verticaux sont les niveaux du modele
39 !           ZCDL si les niveaux verticaux sont des niveaux Z=constante donnes au programme
40 !           PCDL si les niveaux verticaux sont des niveaux P=constante donnes au programme
41 !
42 !  pour les formats *CDL,*Z*,*P*, 2 types de grille horizontale sont possibles:
43 !    'CONF' grille reguliere sur le plan de projection (conforme ou cartesien)
44 !    'LALO' grille reguliere en lat-lon
45 !             dans ce cas les composantes du vent sont transformees
46 !             en composantes zonales et méridiennes.
47 !!
48 !!    EXTERNAL
49 !!    --------
50 !!          FROM_COMPUTING_UNITS: retour aux unites initiales  avant ecriture
51 !!                               = passage inverse a celui realise par
52 !!                                 TO_COMPUTING_UNITS      
53 !!          appele par writevar,writecdl,writellhv 
54 !!              et par extractdia avant l ecriture au format FREE
55 !!    REFERENCE
56 !!    ---------
57 !!
58 !!    AUTHORS
59 !!    -------
60 !!    I. Mallet , N. Asencio, J. Stein
61 !!
62 !!
63 !!    MODIFICATIONS
64 !!    -------------
65 !!      Original    17/03/2003
66 !       call to dd and ff routines
67 !       call to writeLLHV if LLHV
68 !       clean writevar to delete choice LLHV inside this routine
69 !       add PCDL,LLZV,llzv,LLPV,llpv cases
70 !       allow a zoom 0,0,jdeb,jfin or ideb,ifin,0,0 or 0,0,0,0  05/2005
71 !        add ALT 3Dfield if KCDL, add the LAT and LON 3Dfields if CONF and *CDL
72 !-------------------------------------------------------------------------------
73 !
74 !*       0.    DECLARATIONS
75 !              ------------
76 !
77 ! modules MesoNH
78 USE MODD_CONF, ONLY: NVERB
79 USE MODD_PARAMETERS, ONLY: JPHEXT,JPVEXT,XUNDEF
80 USE MODD_DIM1, ONLY: NIMAX,NJMAX,NKMAX
81 USE MODD_GRID, ONLY: XLATORI,XLONORI
82 USE MODD_GRID1, ONLY: XZS,XZZ,XLAT,XLON,XXHAT,XYHAT
83 USE MODD_LUNIT1, ONLY: CLUOUT
84 USE MODE_GRIDPROJ  ! subroutines SM_XYHAT et SM_LATLON 
85 USE MODI_UV_TO_ZONAL_AND_MERID
86 USE MODI_HOR_INTERP_4PTS
87 USE MODI_ZINTER
88 USE MODI_PINTER                          
89 ! modules DIACHRO
90 USE MODD_FILES_DIACHRO
91 USE MODN_NCAR,  ONLY: XSPVAL  
92 USE MODD_ALLOC_FORDIACHRO, ONLY: XVAR, &         ! XVAR(i,j,k,t,n,p)
93                                  XTRAJZ, &       ! XTRAJZ(k,t,n)
94                                  XDATIME, &      ! XDATIME(16,t)
95                                  CTITRE, CUNITE,&! CTITRE(p),CUNITE(p)
96 !* UPG irina
97                                  XTRAJT, &       ! XTRAJT(t,n)
98 !* UPG irina
99                                  NGRIDIA, & ! NGRIDIA(p)
100                                  NGRID
101 USE MODD_COORD, ONLY: XXX,XXY,XXZS, & !  XXX(:,1:7), XXY(:,1:7), XXZS(:,:,1:7)
102                       XXDXHAT,XXDYHAT ! XXDXHAT(:,1:7), XXDYHAT(:,1:7)
103 USE MODD_RESOLVCAR, ONLY: CGROUP, NVERBIA, &
104                           NNDIA, NPROCDIA, NBPROCDIA !pour appel a interp_grids
105 USE MODD_TYPE_AND_LH, ONLY: NIL,NIH,NJL,NJH,NKL,NKH,CTYPE,LICP,LJCP
106 ! modules tools
107 USE MODI_CHANGE_A_GRID 
108 USE MODI_LOW2UP 
109 USE MODI_CREATLINK 
110 USE MODI_DD
111 USE MODI_FF
112 USE MODI_WRITEDIR                                 
113 USE MODI_WRITELLHV                                 
114 USE MODI_WRITECDL                                 
115 USE MODI_WRITEVAR                                 
116 USE MODI_FROM_COMPUTING_UNITS
117 USE MODD_READLH
118 !                                 
119 IMPLICIT NONE                       
120 !
121 !*       0.1   Local variables declarations
122 !
123 INTEGER           :: I
124 INTEGER           :: ILUDIR,IRESP
125 INTEGER           :: JLOOP,JI,JJ,JK,J5,J6,J4,JA,JGR
126 ! zoom lu pour les 6 dimensions possibles
127 INTEGER           :: iideb,iifin,ijdeb,ijfin,ikdeb,ikfin
128 REAL              :: zideb,zifin,zjdeb,zjfin
129 INTEGER, dimension(2) :: iloc
130 INTEGER           :: itinf,itsup,itrajinf,itrajsup,iprocinf,iprocsup
131 ! zoom recalcule en fonction des dimensions du champ traite
132 INTEGER           :: ivarideb,ivarifin,ivarjdeb,ivarjfin,ivarkdeb,ivarkfin
133 INTEGER           :: ivartinf,ivartsup,ivartrajinf,ivartrajsup,ivarprocinf,ivarprocsup
134 INTEGER           :: ivarzmin,ivarzmax
135 INTEGER           :: inbvertz,IND_VERT,IND_LL
136 REAL , allocatable, dimension(:,:,:):: ZWORK3D,ZWORK3D2,zffvent,zdirvent
137 REAL , allocatable, dimension(:,:)  :: zwork2d,zwork2d2
138 REAL , allocatable, dimension(:,:)  :: ZLAT,ZLON
139 ! pour traiter les champs budget deja zoomes
140 REAL , allocatable, dimension(:,:,:,:,:,:):: ZVARSAVE                                 
141 ! pour l interpolation verticale a P=cst : pinter
142 REAL , allocatable, dimension(:,:,:) :: ZPABS                          
143 ! pour les interpolations verticales a P ou Z=cst
144 REAL , allocatable, dimension(:,:,:) :: ZVARZCST
145 REAL , allocatable, dimension(:) :: zlistevert
146 INTEGER :: ikdebzint ! premier niveau a traiter
147 !  pour l interpolation sur grille reguliere lat lon
148 REAL , allocatable, dimension(:,:) :: ZNEWLAT,ZNEWLON,ZNEWX,ZNEWY
149 REAL              :: ZDELTALAT,ZDELTALON
150 !* UPG irina
151 REAL              :: ZLAG
152 !* UPG irina
153 REAL :: zmini,zmaxi
154 INTEGER           :: inetadd ! compteur de champs supp dans le fichier Netcdf
155 INTEGER           :: IFLAGzcst,IGRID
156 INTEGER           :: IDIM1,IDIM2,I1,I2,IZOOMIDEB,IZOOMIFIN,IZOOMJDEB,IZOOMJFIN
157 INTEGER           :: IAN,IMOIS,IJOUR,IHEURE,IMINUTE,ISECONDE
158
159 INTEGER           :: ilocverbia,iret,iret2,iskip,ISAVENGRIDIA,iarg,INDX,IK
160 CHARACTER(LEN=3)  :: YK
161 !                    flag pour initialiser/ne pas initialiser le zoom d
162 !                      d ecriture : 
163 !                      ne pas initialiser quand ajout par le programme
164 !                      des champs ALT LAT LON qui doivent conserver le
165 !                      zoom de l utilisateur
166 INTEGER           :: ino_init_zoom
167 ! **** la taille des variables caracteres contenant les noms
168 !      de fichiers est obligatoirement de 28 ****
169 CHARACTER(LEN=28) :: YFILEIN,YFILEOUT
170 ! **** la longueur du nom ne doit pas depasser 13 car. si le fichier
171 ! contient des groupes a un seul PROCessus, ou 9 si plusieurs PROCessus ****
172 CHARACTER(LEN=13) :: YGROUP,YGROUP_OLD
173 CHARACTER(LEN=20) :: YGROUP_SAVE
174 CHARACTER(LEN=4)  :: YTYPEOUT
175 CHARACTER(LEN=1)  :: YTYPEOUT3
176 CHARACTER(LEN=3)  :: YSUFFIX_file
177 CHARACTER(LEN=250):: YFMTFREE   ! format ecriture des champs si YTYPEOUT='FREE'
178 CHARACTER(LEN=45) :: YFILEOUTFREE ! nom du fichier de sortie si YTYPEOUT='FREE'
179
180 CHARACTER(LEN=5)  :: YFLAGREADVAR ,YFLAGWRITE
181 CHARACTER(LEN=4)  :: YOUTGRID  ! grille en sortie:
182                   !CONF pour rester dans le plan conforme,
183                   ! (le logiciel graphique devra réaliser la projection)
184                   !LALO pour passer à lat,lon réguliers
185 CHARACTER(LEN=28) :: YDUMMYFILE
186 CHARACTER(LEN=11) :: YLUDIR      !  Name of the dir file
187 REAL   , DIMENSION(:,:)  ,ALLOCATABLE        :: ZX,ZY                  
188 !-------------------------------------------------------------------------------
189 !
190 !*       1.     INIT
191 !               ----
192 inetadd=0  !compteur de champs supp dans le fichier Netcdf
193 !
194 !Prints : 0=mini 1=debug mode in extractdia, readvar and writevar , writecdl, writellhv
195 !                3=debug mode in routines diachro'
196 ! nverbia= controle des prints dans les routines diachro
197 ilocverbia=0
198
199 ! dans mesonh Xundef est utilise =999.
200 ! dans les routines diachro XSPVAL est utilisé                  
201 XSPVAL=XUNDEF                                    
202 !
203 ! ouverture d un fichier dir ou vont s ecrire les entrees clavier
204 YLUDIR='dirextract'
205 CALL FMATTR(YLUDIR,YLUDIR,ILUDIR,IRESP)
206 OPEN(UNIT=ILUDIR,FILE=YLUDIR,FORM='FORMATTED')
207 !
208 ! Possibilite de definir un zoom d ecriture 
209 !  definition locale du zoom pour extractdia et writevar, writecdl, writellhv
210 iideb=0
211 iifin=0
212 ijdeb=0
213 ijfin=0
214 ikdeb=0
215 ikfin=0
216 itinf=0
217 itsup=0
218 itrajinf=0
219 itrajsup=0
220 iprocinf=0
221 iprocsup=0
222 !
223 !-------------------------------------------------------------------------------
224 !
225 !*       2.     INPUT FILE AND FORMAT
226 !               ---------------------
227 !
228 !*       2.1   name of file and output format
229 !              ------------------------------
230 !
231 PRINT*, '- Name of the diachro file (without .lfi) ?'
232 READ(5,'(A28)') YFILEIN
233 CALL WRITEDIR(ILUDIR,YFILEIN)
234 !
235 PRINT*, '- type of the output file (DIAC/llhv/llzv/llpv/LLHV/LLZV/LLPV/FREE/KCDL/ZCDL/PCDL)'
236 READ(5,'(A4)')YTYPEOUT
237 CALL WRITEDIR(ILUDIR,YTYPEOUT)
238 PRINT*,'the file ',TRIM(YFILEIN),' will be converted in type ',YTYPEOUT
239 !
240 PRINT*, '- Prints : 0=mini 1=debug mode in extractdia'
241 PRINT*, '                  3=debug mode in routines diachro'
242 PRINT*, '?'
243 READ(5,*)ilocverbia
244 CALL WRITEDIR(ILUDIR,ilocverbia)
245 PRINT*, ' output prints= ',ilocverbia
246 if ( ilocverbia > 2) nverbia=ilocverbia   ! verbosity of diachro routines
247 NVERB=ilocverbia                          ! verbosity of mesonh routines
248 !
249 !*       2.2   other parameters
250 !              ----------------
251 !
252 SELECT CASE (YTYPEOUT)                                   
253   CASE('LLHV','llhv','DIAC','FREE','KCDL','ZCDL','PCDL','llzv','LLZV','llpv','LLPV') ! lecture des choix de l utilisateur
254 !* UPG irina
255       IF ( YTYPEOUT == 'DIAC' ) THEN
256         PRINT*, 'valeur temporelle a ajouter a XTRAJT ? '
257         read(5,*) ZLAG
258         print*,ZLAG
259       ENDIF
260 !* UPG irina
261     IF ( YTYPEOUT == 'FREE' ) THEN
262       PRINT*, '- format of writing for fields ? '
263       PRINT*, '    (fortran syntaxe of FMT in WRITE)'
264       PRINT*,'exemple: (10F9.3) or (8F0.3)'
265       PRINT*, '?'
266       READ(5,'(A)') YFMTFREE
267       CALL WRITEDIR(ILUDIR,YFMTFREE)
268       PRINT*, ' format=', TRIM(YFMTFREE)
269     ENDIF
270     ! lecture du zoom
271     IND_VERT= INDEX(YTYPEOUT(1:4),'Z') + INDEX(YTYPEOUT(1:4),'P') + &
272               INDEX(YTYPEOUT(1:4),'z') + INDEX(YTYPEOUT(1:4),'p')
273     IND_LL= INDEX(YTYPEOUT(1:2),'L') + INDEX(YTYPEOUT(1:2),'l') 
274     IF (IND_LL==0) THEN
275       IF (IND_VERT/=0) THEN
276         ! cas 'ZCDL','PCDL'
277         PRINT*, '- zoom on the 2 first dimensions: '
278         PRINT*, '              ideb,ifin,jdeb,jfin'
279         PRINT*, '0,0,0,0 for the whole physical domain'
280         PRINT*, '-1,-1,-1,-1 for the whole domain'
281         PRINT*, '?'
282         READ(5,*) iideb,iifin,ijdeb,ijfin
283         CALL WRITEDIR(ILUDIR,iideb)
284         CALL WRITEDIR(ILUDIR,iifin)
285         CALL WRITEDIR(ILUDIR,ijdeb)
286         CALL WRITEDIR(ILUDIR,ijfin)
287       ELSE 
288         ! cas 'DIAC','FREE','KCDL'
289         PRINT*, '- zoom on the 3 first dimensions: '
290         PRINT*, '              ideb,ifin,jdeb,jfin,kdeb,kfin'
291         PRINT*, '0,0,0,0,0,0 for the whole physical domain'
292         PRINT*, '-1,-1,-1,-1,-1,-1 for the whole domain'
293         PRINT*, '?'
294         READ(5,*) iideb,iifin,ijdeb,ijfin,ikdeb,ikfin
295         CALL WRITEDIR(ILUDIR,iideb)
296         CALL WRITEDIR(ILUDIR,iifin)
297         CALL WRITEDIR(ILUDIR,ijdeb)
298         CALL WRITEDIR(ILUDIR,ijfin)
299         CALL WRITEDIR(ILUDIR,ikdeb)
300         CALL WRITEDIR(ILUDIR,ikfin)
301       END IF
302     ELSE
303       ! cas 'llzv','LLZV','llpv','LLPV','llhv','LLHV'
304       PRINT*, '- zoom on the 2 first directions: '
305       PRINT*, '              lonmin,lonmax,latmin,latmax'
306       PRINT*, '0.,0.,0.,0. for the whole physical domain'
307       PRINT*, '-1.,-1.,-1.,-1. for the whole domain'
308       PRINT*, '?'
309       READ(5,*) zideb,zifin,zjdeb,zjfin
310       CALL WRITEDIR(ILUDIR,zideb)
311       CALL WRITEDIR(ILUDIR,zifin)
312       CALL WRITEDIR(ILUDIR,zjdeb)
313       CALL WRITEDIR(ILUDIR,zjfin)
314       if(zideb==0. .AND. zifin==0.) then
315         iideb=0 ; iifin=0
316       else if(zideb==-1. .AND. zifin==-1.) then
317         iideb=-1 ; iifin=-1
318       else
319         iideb=-2 ; iifin=-2
320       endif
321       if(zjdeb==0. .AND. zjfin==0.) then
322         ijdeb=0 ; ijfin=0
323       else if(zjdeb==-1. .AND. zjfin==-1.) then
324         ijdeb=-1 ; ijfin=-1
325       else
326         ijdeb=-2 ; ijfin=-2
327       endif
328       !! O.Nuissier
329       !!iideb=zideb ; iifin=zifin ; ijdeb=zjdeb ; ijfin=zjfin
330       !! O.Nuissier
331       IF (IND_VERT==0) THEN
332         ! cas 'llhv','LLHV'
333         PRINT*, '- zoom on the 3rd dimension: '
334         PRINT*, '                 kdeb,kfin'
335         PRINT*, '0,0 for the whole physical domain'
336         PRINT*, '-1,-1 for the whole domain'
337         PRINT*, '?'
338         READ(5,*) ikdeb,ikfin
339         CALL WRITEDIR(ILUDIR,ikdeb)
340         CALL WRITEDIR(ILUDIR,ikfin)
341       END IF
342     END IF
343     PRINT*, '- zoom on the 3 last dimensions : '
344     PRINT*, '   itinf,itsup,itrajinf,itrajsup,iprocinf,iprocsup'
345     PRINT*, '0,0,0,0,0,0 for the whole last dimensions'
346     PRINT*, '?'
347     READ(5,*) itinf,itsup,itrajinf,itrajsup,iprocinf,iprocsup
348     CALL WRITEDIR(ILUDIR,itinf)
349     CALL WRITEDIR(ILUDIR,itsup)
350     CALL WRITEDIR(ILUDIR,itrajinf)
351     CALL WRITEDIR(ILUDIR,itrajsup)
352     CALL WRITEDIR(ILUDIR,iprocinf)
353     CALL WRITEDIR(ILUDIR,iprocsup)
354     IF ((iideb==-2) .AND. (ijdeb==-2)) THEN
355       PRINT*, ' zoom= ',zideb,zifin,zjdeb,zjfin,ikdeb,ikfin&
356                       ,itinf,itsup,itrajinf,itrajsup,iprocinf,iprocsup
357     ELSE
358       PRINT*, ' zoom= ',iideb,iifin,ijdeb,ijfin,ikdeb,ikfin&
359                       ,itinf,itsup,itrajinf,itrajsup,iprocinf,iprocsup
360     END IF
361     IF (IND_VERT/=0) THEN
362       PRINT*, '- Number of vertical levels for ',YTYPEOUT(IND_VERT:IND_VERT),' interpolation ?'
363       READ(5,*) inbvertz
364       CALL WRITEDIR(ILUDIR,inbvertz)
365       PRINT*, '- List of these levels (in meters or in hPa): exemple 500 1500 ?'
366       allocate (zlistevert(inbvertz))
367       READ(5,*) zlistevert
368       DO JI=1,inbvertz
369         CALL WRITEDIR(ILUDIR,zlistevert(JI))
370       END DO
371       PRINT*, ' interpolation for the following ',YTYPEOUT(IND_VERT:IND_VERT),' levels='
372       PRINT*, zlistevert 
373     ENDIF
374     YOUTGRID='CONF'
375     IF (YTYPEOUT/='DIAC' .AND. YTYPEOUT/='llhv' .AND. YTYPEOUT/='LLHV') THEN
376       PRINT *,'- Fields in regular LAt/LOn grid'
377       PRINT *,'  or    in regular grid on CONFormal plan (native MesoNH grid) ?'
378       PRINT *,'LALO/CONF ?'
379       READ(5,*) YOUTGRID
380       CALL WRITEDIR(ILUDIR,YOUTGRID)
381       PRINT*, ' Output grid= ', YOUTGRID
382       PRINT*, ''
383       YSUFFIX_file=YTYPEOUT(1:2)//YTYPEOUT(4:4)
384       IF ( YTYPEOUT(2:4) == 'CDL') THEN
385         PRINT*, '!!!!!!!! Warning !!!!!!!!'
386         PRINT*, 'For the CDL type, the dimensions are initialised'
387         PRINT*, ' with those of the first field:'
388         PRINT*, 'the values of the 6 dimensions must be the maximum that'
389         PRINT*, ' will be treated '
390         PRINT*, '!!!!!!!! Warning !!!!!!!!'
391         PRINT*, 'For the CDL type, the coordinates must be the same'
392         PRINT*, ' for all fields'
393         PRINT*, '(stored in the output file with LAT/LON/ALT groups)'
394         PRINT*, '!!!!!!!!'
395       ENDIF
396     ENDIF
397   CASE DEFAULT
398     PRINT*, 'Incorrect value for the output type:',YTYPEOUT
399     PRINT*, ' the following ones are currently available : DIAC,LLHV,llhv,FREE,KCDL,ZCDL,PCDL,llzv,LLZV,llpv,LLPV'
400     STOP
401 END SELECT
402
403 !*       2.3   init for input file and output file
404 !              -----------------------------------
405 ! in READVAR, input file must be opened before reading
406 YFLAGREADVAR='OPE'
407 ! in WRITE routine, output file is new
408 YFLAGWRITE='NEW'
409
410 !*       2.4   lecture de la pression pour interpolation
411 !              -----------------------------------------
412 IF (INDEX(YTYPEOUT(1:4),'p')/=0 .OR. INDEX(YTYPEOUT(1:4),'P')/=0 )THEN
413   CALL READVAR('PABSM',YFILEIN,YFLAGREADVAR,ilocverbia,iret)
414   IF ( iret /= 0 ) then
415     print *, '- PABSM not found, name of the pressure variable ? '
416     read *,YGROUP
417     CALL WRITEDIR(ILUDIR,YGROUP)
418     CALL READVAR(YGROUP,YFILEIN,YFLAGREADVAR,ilocverbia,iret)
419     IF ( iret /= 0 ) then
420       print *,' interpolation at P=cst not possible because PABSM and ',TRIM(YGROUP),' are not available'
421       STOP
422     ENDIF
423   ENDIF
424   ! stockage de ZPABS utilise par pinter
425   ALLOCATE ( ZPABS(SIZE(XVAR,1),SIZE(XVAR,2),SIZE(XVAR,3))) 
426   ZPABS(:,:,:)=XVAR(:,:,:,1,1,1) 
427 ENDIF
428 !
429 !-------------------------------------------------------------------------------
430 !
431 !*       3.    LOOP ON GROUPS IN THE FILE
432 !              --------------------------
433 !
434 DO JGR=1,10000
435   !    
436   !*      3.0  preparation pour la lecture du champ suivant
437   !
438   ino_init_zoom=0
439   PRINT*,'- Name of the group in upper case (13 characters max.)'
440   PRINT*,' (ex: THM or DD or FF or DD10 or FF10 or LAT or LON or ALT)'
441   PRINT*,'(GROUP for the list of groups, END to stop)?'
442   READ(5,'(A13)',END=88) CGROUP
443   CALL WRITEDIR(ILUDIR,CGROUP)
444   CGROUP=ADJUSTL(CGROUP)
445   CALL LOW2UP(CGROUP)
446   IF (CGROUP=='END') GO TO 88
447   ! point de reprise pour forcer l ecriture des champs ALT,LAT,LON 
448   ! dans les fichiers netcdf
449 77 CONTINUE
450   YGROUP_SAVE=CGROUP(1:13)
451   YK=''
452   INDX=INDEX(CGROUP,'_K_')
453   IF (INDX/=0) THEN
454     CGROUP=YGROUP_SAVE(1:INDX-1)
455     YK(1:3)=YGROUP_SAVE(INDX+3:INDX+5)
456     READ(YK,'(I3)') IK
457   END IF
458   IF (CGROUP(1:5)/='GROUP') &
459     PRINT*,'you asked for the following record: ',TRIM(CGROUP)
460   !
461   !*      3.1  Lecture et initialisation du tableau XVAR
462   !            passé en module MODD_ALLOC_FORDIACHRO
463   !
464   !
465   !      3.1.1 Cas particulier pour le vent
466   !
467   IF ( CGROUP(1:2) == 'UM' .OR. &
468        CGROUP(1:2) == 'VM' .OR. &
469        CGROUP(1:2) == 'DD' .OR. &
470        CGROUP(1:2) == 'FF'      )  THEN
471     !
472     IF ( (CGROUP(1:2)=='UM'.OR.CGROUP(1:2)=='VM') .AND. &
473           YOUTGRID(1:4) /= 'LALO'                       ) THEN
474       ! Lecture du champ U ou V sans calcul 
475       ! les composantes du vent restent dans le plan conforme
476       CALL READVAR(CGROUP,YFILEIN,YFLAGREADVAR,ilocverbia,iret)
477     ELSE
478       ! Lecture des 2 composantes du vent 
479       !(stockees dans les tableaux ZWORK3D et ZWORK3D2)
480       IF (LEN(TRIM(CGROUP)) ==2) THEN
481         YGROUP='UM'
482       ELSE IF (LEN(TRIM(CGROUP)) ==4) THEN
483         YGROUP='UM'//CGROUP(3:4)
484       ELSE
485         print*,'** problem with the name of group: ',CGROUP
486         CYCLE
487       ENDIF
488       CALL READVAR(YGROUP,YFILEIN,YFLAGREADVAR,ilocverbia,iret)
489       IF ( iret /= 0 ) then
490         print *,TRIM(CGROUP),': ',TRIM(YGROUP),' not available'
491         IF (LEN(TRIM(CGROUP)) ==2) THEN
492           YGROUP='UT'
493         ELSE IF (LEN(TRIM(CGROUP)) ==4) THEN
494           YGROUP='UT'//CGROUP(3:4)
495         ENDIF
496         CALL READVAR(YGROUP,YFILEIN,YFLAGREADVAR,ilocverbia,iret2)
497         IF ( iret2 /= 0 ) then
498           print *,'** no processing for ',TRIM(CGROUP), &
499                   ' because UM and ',TRIM(YGROUP),' are not available'
500           CYCLE
501         ENDIF
502       ENDIF
503       ! allocation du tableau de stockage de la 1e composante du vent
504       ALLOCATE(ZVARSAVE(size(XVAR,1),size(XVAR,2),size(XVAR,3),   &
505                         size(XVAR,4),size(XVAR,5),size(XVAR,6))   )
506       ZVARSAVE=XVAR
507       !
508       IF (LEN(TRIM(CGROUP)) ==2) THEN
509         YGROUP='VM'
510       ELSE IF (LEN(TRIM(CGROUP)) ==4) THEN
511         YGROUP='VM'//CGROUP(3:4)
512       ENDIF
513       CALL READVAR(YGROUP,YFILEIN,YFLAGREADVAR,ilocverbia,iret)
514       IF ( iret /= 0 ) then
515         print *,TRIM(CGROUP),': ',TRIM(YGROUP),' not available'
516         IF (LEN(TRIM(CGROUP)) ==2) THEN
517           YGROUP='VT'
518         ELSE IF (LEN(TRIM(CGROUP)) ==4) THEN
519           YGROUP='VT'//CGROUP(3:4)
520         ENDIF
521         CALL READVAR(YGROUP,YFILEIN,YFLAGREADVAR,ilocverbia,iret2)
522         IF ( iret2 /= 0 ) then
523           print *,'** no processing for ',TRIM(CGROUP), &
524                   ' because VM and ',TRIM(YGROUP),' are not available'
525           CYCLE
526         ENDIF
527         iret=iret2
528       ENDIF
529       !
530       ! Calcul de ff
531       IF (CGROUP(1:2) == 'FF' ) THEN
532         IF (LEN(TRIM(CGROUP)) ==2) THEN
533           YGROUP='VENTFF'
534         ELSE IF (LEN(TRIM(CGROUP)) ==4) THEN
535           YGROUP='VENT'//CGROUP(3:4)//'FF'
536         ENDIF
537         ! allocation du tableau de calcul
538         IF (allocated(ZWORK3D)) DEALLOCATE(ZWORK3D)
539         ALLOCATE(ZWORK3D(size(XVAR,1),size(XVAR,2),size(XVAR,3)))
540         ZWORK3D(:,:,:)=XSPVAL
541         DO J6=1,SIZE(XVAR,6)
542           IGRID=NGRIDIA(J6)
543           DO J5=1,SIZE(XVAR,5)
544           DO J4=1,SIZE(XVAR,4)
545             CALL FF (ZVARSAVE(:,:,:,J4,J5,J6),XVAR(:,:,:,J4,J5,J6),ZWORK3D, &
546                      JPVEXT,JPHEXT,IGRID)
547             XVAR(:,:,:,J4,J5,J6)=ZWORK3D(:,:,:)
548           END DO
549           END DO
550           ! initialisation des variables necessaires a l ecriture
551           CGROUP=YGROUP
552           CTITRE(J6)=YGROUP
553           NGRIDIA(J6)=1
554         END DO
555         DEALLOCATE(ZWORK3D)
556         ! Calcul de dd par rapport au Nord geographique
557       ELSE IF (CGROUP(1:2) == 'DD') THEN
558         IF (CTYPE=='CART' .OR. CTYPE=='MASK' .OR. CTYPE=='SPXY') THEN 
559           IF (LEN(TRIM(CGROUP)) ==2) THEN
560             YGROUP='VENTDD'
561           ELSE IF (LEN(TRIM(CGROUP)) ==4) THEN
562             YGROUP='VENT'//CGROUP(3:4)//'DD'
563           ENDIF
564           ! allocation du tableau de calcul
565           IF (allocated(ZWORK3D)) DEALLOCATE(ZWORK3D)
566           ALLOCATE(ZWORK3D(size(XVAR,1),size(XVAR,2),size(XVAR,3)))
567           DO J6=1,SIZE(XVAR,6)
568             IGRID=NGRIDIA(J6)
569             DO J5=1,SIZE(XVAR,5)
570             DO J4=1,SIZE(XVAR,4)
571               iskip=1 ! tous les points de grille
572               CALL DD(ZVARSAVE(:,:,:,J4,J5,J6),XVAR(:,:,:,J4,J5,J6),ZWORK3D, &
573                       iskip,IGRID,PLON=XLON(NIL:NIH,NJL:NJH))
574               XVAR(:,:,:,J4,J5,J6)=ZWORK3D(:,:,:)
575             END DO
576             END DO
577             ! initialisation des variables necessaires a l ecriture
578             CGROUP=YGROUP
579             CTITRE(J6)=YGROUP
580             CUNITE(J6)='degrees'
581             NGRIDIA(J6)=1
582           END DO
583           DEALLOCATE(ZWORK3D)
584         ELSE
585           print *,'** processing of ',TRIM(CGROUP),' is not performed for CTYPE= ',CTYPE
586           CYCLE
587         ENDIF
588       ELSE IF (CGROUP(1:2) == 'UM' .OR. CGROUP(1:2) == 'VM') THEN
589         IF (CTYPE=='CART' .OR. CTYPE=='MASK' .OR. CTYPE=='SPXY') THEN 
590         ! Calcul des composantes zonale et meridienne
591         !(YOUTGRID(1:4) == 'LALO') avec la routine UV_TO_ZONAL_AND_MERID
592           print*,' Translate to meridional and zonal wind components'
593           ALLOCATE(ZWORK3D(size(XVAR,1),size(XVAR,2),size(XVAR,3)))
594           ALLOCATE(ZWORK3D2(size(XVAR,1),size(XVAR,2),size(XVAR,3)))
595           IF (ilocverbia >= 3 ) then
596             print *,'before UV_TO_ZONAL_AND_MERID KGRID=23'
597             print *,' dimensions of the input arrays',size(ZVARSAVE,1),&
598                                       size(ZVARSAVE,2),size(ZVARSAVE,3)
599             print *,size(XVAR,1),size(XVAR,2),size(XVAR,3)
600             print *,' dimensions of the output arrays',size(ZWORK3D,1),&
601                                        size(ZWORK3D,2),size(ZWORK3D,3)
602             print *,size(ZWORK3D2,1),size(ZWORK3D2,2),size(ZWORK3D2,3)
603           ENDIF
604           DO J6=1,SIZE(XVAR,6)
605             IGRID=NGRIDIA(J6)
606             DO J5=1,SIZE(XVAR,5)
607             DO J4=1,SIZE(XVAR,4)
608               CALL UV_TO_ZONAL_AND_MERID(ZVARSAVE(:,:,:,J4,J5,J6), &
609                                          XVAR(:,:,:,J4,J5,J6),     &
610                                          23,PZC=ZWORK3D,PMC=ZWORK3D2)
611               IF (CGROUP(1:1) == 'U' ) THEN
612                 XVAR(:,:,:,J4,J5,J6)=ZWORK3D(:,:,:)
613               ENDIF
614               IF (CGROUP(1:1) == 'V' ) THEN
615                 XVAR(:,:,:,J4,J5,J6)=ZWORK3D2(:,:,:)
616               ENDIF
617             END DO
618             END DO
619           END DO
620           IF (ilocverbia >= 3 ) then
621             print *,'after UV_TO_ZONAL_AND_MERID KGRID=23'
622           END IF
623           ! Stockage dans le tableau XVAR qui est le tableau ecrit
624           ! de la composante souhaitée
625           IF (CGROUP(1:1) == 'U' ) THEN
626             print *, ' U zonal wind component'
627             IF (LEN(TRIM(CGROUP)) ==2) THEN
628               YGROUP='UZON'
629             ELSE IF (LEN(TRIM(CGROUP)) ==4) THEN
630               YGROUP='U'//CGROUP(3:4)//'ZON'
631             ENDIF
632             CGROUP=YGROUP
633             CTITRE(:)='U zonal wind component'
634           ENDIF
635           IF (CGROUP(1:1) == 'V' ) THEN
636             print *, ' V meridian wind component'
637             IF (LEN(TRIM(CGROUP)) ==2) THEN
638               YGROUP='VMED'
639             ELSE IF (LEN(TRIM(CGROUP)) ==4) THEN
640               YGROUP='V'//CGROUP(3:4)//'MED'
641             END IF
642             CGROUP=YGROUP
643             CTITRE(:)='V meridian wind component'
644           ENDIF
645           DEALLOCATE(ZWORK3D,ZWORK3D2)
646         ELSE
647           print *,' No processing of UZON and VMED for CTYPE= ',CTYPE
648           CYCLE
649         ENDIF
650       ENDIF
651       DEALLOCATE(ZVARSAVE)
652     ENDIF
653   !
654   !      3.1.2 LATitude ou LONgitude de chaque point de la grille conforme
655   !
656   ELSE IF (CGROUP(1:3)=='LAT' .OR. CGROUP(1:3)=='LON') THEN
657     print *, 'LAT/LON asked and YFLAGREADVAR=', YFLAGREADVAR
658    IF ( YFLAGREADVAR /= 'NOP') THEN
659     ! Lecture d un champ 2D quelconque pour initialiser XLAT et XLON
660     CALL READVAR('ZSBIS',YFILEIN,YFLAGREADVAR,ilocverbia,iret)
661     IF ( iret /= 0 ) then
662      ! cas de fichier diachronique sans ZSBIS
663       print *, '- Name of one group in upper case '
664       read *,YGROUP
665       CALL WRITEDIR(ILUDIR,YGROUP)
666       CALL LOW2UP(YGROUP)
667       CALL READVAR(YGROUP,YFILEIN,YFLAGREADVAR,ilocverbia,iret)
668       IF ( iret /= 0 ) then
669          print * ,'**group ', TRIM(YGROUP) , 'not found'
670          stop
671       ENDIF
672     ENDIF
673    ENDIF
674     ! init du tableau XVAR au champ souhaite
675     DEALLOCATE(XVAR)
676     ALLOCATE(XVAR(size(XLAT,1),size(XLAT,2),1,1,1,1) )
677     IF (CGROUP(1:3)=='LAT') THEN
678       XVAR(:,:,1,1,1,1)=XLAT(:,:)
679       CTITRE(1)='latitudes'
680       CUNITE(1)='degrees_north'
681     ELSE IF (CGROUP(1:3)=='LON') THEN
682       XVAR(:,:,1,1,1,1)=XLON(:,:)
683       CTITRE(1)='longitudes'
684       CUNITE(1)='degrees_east'
685     ENDIF
686   !
687   !      3.1.3 ALTitude de chaque point de la grille conforme
688   !
689   ELSE IF (CGROUP(1:3)=='ALT') THEN
690     print *, 'ALT asked and YFLAGREADVAR=', YFLAGREADVAR
691     IF(CTYPE=='SSOL'.OR.CTYPE=='DRST'.OR.CTYPE=='RAPL'.OR.CTYPE=='RSPL') THEN
692       IF ( YFLAGREADVAR == 'NOP') THEN
693       ! altitude des niveaux du groupe precedent dans XTRAJZ
694         print *,'warning, for CTYPE=',CTYPE,' ALTitude of previous group (',TRIM(YGROUP_OLD),')'
695         DEALLOCATE(XVAR)
696         ALLOCATE(XVAR(1,1,size(XTRAJZ,1),1,1,1))
697         XVAR(1,1,:,1,1,1)=XTRAJZ(:,1,1)
698       ELSE
699         print*,'** no processing with ALT at the first group'
700         GOTO 99
701       ENDIF
702     ELSE
703       IF ( YFLAGREADVAR /= 'NOP') THEN
704         ! Lecture d un champ 2D quelconque pour initialiser les tableaux XZZ
705         CALL READVAR('ZSBIS',YFILEIN,YFLAGREADVAR,ilocverbia,iret)
706         IF ( iret /= 0 ) then
707           ! cas de fichier diachronique sans ZSBIS
708           print *, '- Name of one group in upper case '
709           read *,YGROUP
710           CALL WRITEDIR(ILUDIR,YGROUP)
711           CALL LOW2UP(YGROUP)
712           CALL READVAR(YGROUP,YFILEIN,YFLAGREADVAR,ilocverbia,iret)
713           IF ( iret /= 0 ) then
714             print * ,'** group ', TRIM(YGROUP) , 'not found'
715             stop
716           ENDIF
717         ENDIF
718       ENDIF
719       ! init de XZZ a la grille de masse ( par defaut readvar 
720       ! l initialise a la grille 4 des  vitesse verticales W)
721       CALL COMPCOORD_FORDIACHRO(1)
722       ! init du tableau XVAR au champ souhaite
723       DEALLOCATE(XVAR)
724       ALLOCATE(XVAR(size(XZZ,1),size(XZZ,2),size(XZZ,3),1,1,1))
725       XVAR(:,:,:,1,1,1)=XZZ(:,:,:)
726       ! retour au XZZ grille 4
727       CALL COMPCOORD_FORDIACHRO(4)
728     ENDIF
729     CTITRE(1)='model levels altitudes ASL'
730     CUNITE(1)='meters'
731   !
732   !      3.1.4 Default case
733   !
734   ELSE
735     !
736     ! Lecture du  champ CGROUP et stockage dans XVAR
737     ! + Initialisation (si YFLAGREADVAR='OPE') des variables
738     ! des modules (cf USE en debut de programme)
739     ! Appel a menu_diachro pour la liste des groupes si CGROUP(1:5)=='GROUP'
740     !
741     CALL READVAR(CGROUP,YFILEIN,YFLAGREADVAR,ilocverbia,iret)
742     IF (CGROUP(1:5)=='GROUP') CYCLE
743 !* UPG irina
744     print*,'XTRAJT ',SIZE(XTRAJT,1),SIZE(XTRAJT,2)
745     print*,'XTRAJT av. modif. pour les .000 ',XTRAJT(:,:)
746     XTRAJT(:,:)=XTRAJT(:,:)+ZLAG
747 !* UPG irina
748     !
749   ENDIF
750   !
751   IF ( iret == 0 ) THEN
752     zmini=MINVAL(XVAR(:,:,:,:,:,:),MASK=XVAR(:,:,:,:,:,:)/=XSPVAL)
753     zmaxi=MAXVAL(XVAR(:,:,:,:,:,:),MASK=XVAR(:,:,:,:,:,:)/=XSPVAL)
754     print * ,' After read, min,max of the variable ',TRIM(CGROUP),'=', zmini,zmaxi
755     !       
756     !*      3.2  Init des bornes min max du zoom en fonction des
757     !            dimensions du tableau XVAR traite
758     !
759    IF ( ino_init_zoom == 0) THEN
760     IF (iideb == 0 .AND. iifin == 0 ) THEN
761       ivarideb=NREADIL ; ivarifin=NREADIH
762       IF (ivarideb/=ivarifin) THEN  ! domI/=1
763         ivarideb=MAX(1+JPHEXT,NREADIL) 
764         ivarifin=MIN(SIZE(XVAR,1)-JPHEXT,NREADIH)
765         !IF (ivarifin <= 0) THEN
766          ! dimI =1
767         !ivarideb=1 ; ivarifin=SIZE(XVAR,1)
768         !ENDIF
769       ENDIF
770     ELSE IF (iideb == -1 .AND. iifin == -1 ) THEN
771       ivarideb=MAX(1,NREADIL) 
772       ivarifin=MIN(SIZE(XVAR,1),NREADIH)
773     ELSE IF (iideb == -2 .AND. iifin == -2 ) THEN
774       ivarideb=-2
775       iideb=1+JPHEXT
776       IF (zideb >= minval(XLON)) THEN
777         DO JJ=1,SIZE(XLON,2)
778           ivarideb=MAX(MIN(COUNT(XLON(:,JJ)<zideb),SIZE(XLON,1)),iideb)
779           iideb=ivarideb
780         END DO
781       ENDIF
782       ivarifin=-2
783       iifin=1+JPHEXT
784       IF (zifin <= maxval(XLON)) THEN
785         DO JJ=1,SIZE(XLON,2)
786           ivarifin=MAX(MIN(COUNT(XLON(:,JJ)<zifin),SIZE(XLON,1)),iifin)
787           iifin=ivarifin
788         END DO
789       ENDIF
790     ELSE
791       ivarideb=max(iideb,NREADIL)
792       ivarifin=min(iifin,NREADIH)
793       ivarideb=min(ivarideb,ivarifin)
794     ENDIF
795     IF(ijdeb == 0 .AND. ijfin == 0) THEN
796       ivarjdeb=NREADJL ; ivarjfin=NREADJH
797       IF (ivarjdeb/=ivarjfin) THEN  ! domJ/=1
798         ivarjdeb=MAX(1+JPHEXT,NREADJL)
799         ivarjfin=MIN(SIZE(XVAR,2)-JPHEXT,NREADJH)
800         !IF (ivarjfin <= 0) THEN
801          ! dimJ =1
802          !ivarjdeb=1 ; ivarjfin=SIZE(XVAR,2)
803         !ENDIF
804       ENDIF
805     ELSE IF (ijdeb == -1 .AND. ijfin == -1 ) THEN
806       ivarjdeb=MAX(1,NREADJL)
807       ivarjfin=MIN(SIZE(XVAR,2),NREADJH)
808     ELSE IF (ijdeb == -2 .AND. ijfin == -2 ) THEN
809       ivarjdeb=-2
810       ijdeb=1+JPHEXT
811       IF (zjdeb >= minval(XLAT)) THEN
812         DO JI=1,SIZE(XLAT,1)
813           ivarjdeb=MAX(MIN(COUNT(XLAT(JI,:)<zjdeb),SIZE(XLAT,2)),ijdeb)
814           ijdeb=ivarjdeb
815         END DO
816       ENDIF
817       ivarjfin=-2
818       ijfin=1+JPHEXT
819       IF (zjfin <= maxval(XLAT)) THEN
820         DO JI=1,SIZE(XLAT,1)
821           ivarjfin=MAX(MIN(COUNT(XLAT(JI,:)<zjfin),SIZE(XLAT,2)),ijfin)
822           ijfin=ivarjfin
823         END DO
824       ENDIF
825     ELSE
826       ivarjdeb=max(ijdeb,NREADJL)
827       ivarjfin=min(ijfin,NREADJH)
828       ivarjdeb=min(ivarjdeb,ivarjfin)
829     ENDIF
830     IF(ivarideb==-2 .OR. ivarifin==-2 .OR. ivarjdeb==-2 .OR.  ivarjfin==-2) THEN
831       print *,'****zoom provided is not included in the FM-file grid'
832       print *,'LON (zoom: ',zideb,zifin,') (file: ',minval(XLON),maxval(XLON)
833       print *,'LAT (zoom: ',zjdeb,zjfin,') (file: ',minval(XLAT),maxval(XLAT)
834       GOTO 99
835     ENDIF
836     IF (IND_VERT/=0) THEN
837       ivarzmin=1   ; ivarzmax=inbvertz
838     ELSE
839       ivarzmin=MAX(1,NREADKL)  ; ivarzmax=MIN(SIZE(XVAR,3),NREADKH)
840       inbvertz=ivarzmax-ivarzmin+1
841     ENDIF
842     IF (ikdeb == 0 .AND. ikfin == 0 ) THEN
843       ivarkdeb=NREADKL ; ivarkfin=NREADKH
844       IF (ivarkdeb/=ivarkfin) THEN  ! domK/=1
845         ivarkdeb=MAX(1+JPVEXT,NREADKL)
846         ivarkfin=min(ivarzmax,SIZE(XVAR,3)-JPVEXT)
847         !IF (ivarkfin <= 0) THEN
848          ! dimK =1
849          !ivarkdeb=1 ; ivarKfin=SIZE(XVAR,3)
850         !ENDIF
851       ENDIF
852     ELSEIF (ikdeb == -1 .AND. ikfin ==-1 ) THEN
853       ivarkdeb=ivarzmin
854       ivarkfin=ivarzmax
855     ELSE
856       ivarkdeb=max(ikdeb,ivarzmin)
857       ivarkfin=min(ikfin,ivarzmax)
858       ivarkdeb=min(ivarkdeb,ivarkfin)
859     ENDIF   
860     IF (INDX/=0) THEN
861       ivarkdeb=IK ; ivarkfin=IK
862     END IF
863    ENDIF
864
865     IF (itinf == 0 .AND. itsup == 0 ) THEN
866       ivartinf=1 ; ivartsup=SIZE(XVAR,4)
867     ELSE
868       ivartinf=max(itinf,1)
869       ivartsup=min(itsup,SIZE(XVAR,4))
870       ivartinf=min(ivartinf,ivartsup)
871     ENDIF
872     IF (itrajinf == 0 .AND. itrajsup == 0 ) THEN
873       ivartrajinf=1 ; ivartrajsup=SIZE(XVAR,5)
874     ELSE
875       ivartrajinf=max(itrajinf,1)
876       ivartrajsup=min(itrajsup,SIZE(XVAR,5))
877       ivartrajinf=min(ivartrajinf,ivartrajsup)
878     ENDIF
879     IF (iprocinf == 0 .AND. iprocsup == 0 ) THEN
880       ivarprocinf=1 ; ivarprocsup=SIZE(XVAR,6)
881     ELSE
882       ivarprocinf=max(iprocinf,1)
883       ivarprocsup=min(iprocsup,SIZE(XVAR,6))
884       ivarprocinf=min(ivarprocinf,ivarprocsup)
885     ENDIF
886     if (ilocverbia > 0 ) then
887       PRINT*,' Zoom limits initialized with:'
888       PRINT*,'ivarideb,ivarifin,ivarjdeb,ivarjfin,ivarkdeb,ivarkfin',&
889             ivarideb,ivarifin,ivarjdeb,ivarjfin,ivarkdeb,ivarkfin 
890       PRINT*,'ivartinf,ivartsup,ivartrajinf,ivartrajsup,ivarprocinf,ivarprocfin',&
891             ivartinf,ivartsup,ivartrajinf,ivartrajsup,ivarprocinf,ivarprocsup 
892     endif
893     !
894     !*      3.3  Ecriture  du tableau XVAR (module MODD_ALLOC_FORDIACHRO) 
895     !
896     print *,' Write with the format ', YTYPEOUT(1:4)
897     SELECT CASE(YTYPEOUT(1:4))
898       !
899       CASE('DIAC')
900         CALL WRITEVAR(ivarideb,ivarifin,ivarjdeb,ivarjfin,ivarkdeb,ivarkfin,&
901                       ivartinf,ivartsup,ivartrajinf,ivartrajsup,ivarprocinf,ivarprocsup,  &
902                       CGROUP,YFILEIN,YFLAGWRITE,'2  ',ilocverbia,iret)
903         if (ilocverbia > 0 ) then
904           print*,'WRITEVAR return= ',iret
905         end if
906       !
907       CASE('FREE')
908         if (ilocverbia >= 0 ) then
909           print*,' format ',YTYPEOUT
910           print*,' domaine for writting : ideb,ifin,jdeb,jfin,kdeb,kfin', &
911                  ',itinf,itsup,itrajinf,itrajsup,iprocinf,iprocsup= ', &
912               ivarideb,ivarifin,ivarjdeb,ivarjfin,ivarkdeb,ivarkfin,&
913               ivartinf,ivartsup,ivartrajinf,ivartrajsup,ivarprocinf,ivarprocsup 
914         endif
915         !  Retour aux unites initiales si necessaire
916         CALL FROM_COMPUTING_UNITS(CGROUP,CUNITE(1)) 
917         !
918         YFILEOUTFREE=ADJUSTL(ADJUSTR(YFILEIN)//'.'//ADJUSTL(ADJUSTR(CGROUP)))
919         OPEN (UNIT=7,STATUS='NEW',FORM='FORMATTED',FILE=YFILEOUTFREE)
920         ! a. Ecriture de l entete
921         !temps courant
922         IAN=XDATIME(13,1)
923         IMOIS=XDATIME(14,1)
924         IJOUR=XDATIME(15,1)
925         IHEURE=XDATIME(16,1)/3600
926         IMINUTE=(XDATIME(16,1)-(IHEURE*3600))/60
927         ISECONDE=ISECONDE-(IHEURE*3600)-(IMINUTE*60)
928         WRITE(7,*) ivarideb,ivarifin,ivarjdeb,ivarjfin,ivarkdeb,ivarkfin,&
929                    ivartinf,ivartsup,ivartrajinf,ivartrajsup,ivarprocinf,ivarprocsup,&
930                    IAN,IMOIS,IJOUR,IHEURE,IMINUTE ,&
931                    'format ligne1=  12 Indices (.deb .fin) du ',&
932                    'tableau  an mois jour hUTC minute'
933         ! b. ecriture des données au fmt choisi par l utilisateur
934         WRITE(7,FMT=YFMTFREE) &
935          XVAR(ivarideb:ivarifin,ivarjdeb:ivarjfin,ivarkdeb:ivarkfin,&
936               ivartinf:ivartsup,ivartrajinf:ivartrajsup,ivarprocinf:ivarprocsup)
937         PRINT*,'File ',TRIM(YFILEOUTFREE),' available'
938         CLOSE(7)
939       !
940       CASE('LLHV','llhv')
941         CALL WRITELLHV(ivarideb,ivarifin,ivarjdeb,ivarjfin, &
942                        ivarkdeb,ivarkfin,ivartinf,ivartsup, &
943                        ivartrajinf,ivartrajsup,ivarprocinf,ivarprocsup,&
944                        CGROUP,YFILEIN,YFLAGWRITE,YTYPEOUT,&
945                        ilocverbia,iret)       
946         if (ilocverbia > 0 ) then
947           print*,' WRITELLHV return= ',iret
948         end if
949       !
950       CASE('KCDL','ZCDL','PCDL','LLZV','LLPV','llpv','llzv')
951         ! replace field at mass points
952         If (ALLOCATED(ZWORK3D))DEALLOCATE(ZWORK3D)
953         If (ALLOCATED(ZWORK3D2))DEALLOCATE(ZWORK3D2)
954         ALLOCATE(ZWORK3D(SIZE(XVAR,1),SIZE(XVAR,2),SIZE(XVAR,3)))
955         ALLOCATE(ZWORK3D2(SIZE(XVAR,1),SIZE(XVAR,2),SIZE(XVAR,3)))
956         DO J6=ivarprocinf,ivarprocsup
957           IGRID=NGRIDIA(J6)
958           IF(SIZE(XVAR,3)/=1 .OR. IGRID/=4) THEN 
959             ! pas d interpolation verticale pour champ 2D
960             DO J5=ivartrajinf,ivartrajsup
961               DO J4=ivartinf,ivartsup
962                 ZWORK3D(:,:,:)=XVAR(:,:,:,J4,J5,J6)
963                 print *,' mass point grid for J4,J5,J6=',J4,J5,J6
964                 CALL CHANGE_A_GRID(ZWORK3D,IGRID,ZWORK3D2)
965                 ! IGRID=1 en sortie de change_a_grid
966                 XVAR(:,:,:,J4,J5,J6)=ZWORK3D2(:,:,:)
967               ENDDO
968             ENDDO
969           ENDIF
970         ENDDO
971         DEALLOCATE(ZWORK3D,ZWORK3D2)
972         !
973         ! a. reinit avant ecriture de la grille verticale correspondant a la
974         !grille de masse sur laquelle le champ a ete interpole
975         IFLAGzcst=0
976         IF (IND_VERT/=0) THEN
977           IF ( CGROUP == 'ALT' ) THEN
978             ! ecriture de la liste des niveaux verticaux 
979             IFLAGzcst=1
980             DEALLOCATE(XVAR)
981             allocate(XVAR(1,1,inbvertz,1,1,1))
982             XVAR(1,1,:,1,1,1)=zlistevert
983             ivarideb=1 ; ivarifin=1
984             ivarjdeb=1 ; ivarjfin=1
985             ivarkdeb=1 ; ivarkfin=inbvertz
986             CTITRE(1)='vertical_levels'
987             CUNITE(1)='user choice'
988             IF ( YTYPEOUT(IND_VERT:IND_VERT) == 'z' .OR.  YTYPEOUT(IND_VERT:IND_VERT) == 'Z' ) THEN
989               CUNITE(1)='meters'
990             ENDIF
991             IF ( YTYPEOUT(IND_VERT:IND_VERT) == 'p' .OR.  YTYPEOUT(IND_VERT:IND_VERT) == 'P' ) THEN
992               CUNITE(1)='hPa'
993             ENDIF
994           ENDIF
995         ! b. interpolation eventuelle selon la verticale 
996           IF( SIZE(XVAR,3)>1 .AND. SIZE(XVAR,2)>1 .AND. SIZE(XVAR,1)>1 ) THEN
997             ! ALT, LON, LAT et chps 2D ne passent pas cette partie 
998             if (ilocverbia >= 0 ) then
999               print*,' Interpolations on ',inbvertz,' ', &
1000                      YTYPEOUT(IND_VERT:IND_VERT),'-levels'
1001             endif
1002             if (ilocverbia >= 1 .AND. IND_VERT/=0) THEN
1003               print*,'levels= ',zlistevert 
1004             endif
1005             ALLOCATE(ZVARSAVE(size(XVAR,1),size(XVAR,2),size(XVAR,3),   &
1006                               size(XVAR,4),size(XVAR,5),size(XVAR,6))   )
1007             ZVARSAVE=XVAR
1008             ALLOCATE(ZWORK3D(SIZE(XVAR,1),SIZE(XVAR,2),SIZE(XVAR,3)))
1009             ALLOCATE(ZVARZCST(SIZE(XVAR,1),SIZE(XVAR,2),inbvertz))
1010             DEALLOCATE(XVAR)
1011             ALLOCATE(XVAR(SIZE(ZVARSAVE,1),SIZE(ZVARSAVE,2),SIZE(ZVARZCST,3),&
1012                           size(ZVARSAVE,4),size(ZVARSAVE,5),size(ZVARSAVE,6)))
1013             DO J6=ivarprocinf,ivarprocsup
1014               IGRID=NGRIDIA(J6)
1015               ! init du tableau des altitudes  XZZ pour la grille= IGRID
1016               CALL COMPCOORD_FORDIACHRO(IGRID)
1017               DO J5=ivartrajinf,ivartrajsup
1018                 DO J4=ivartinf,ivartsup
1019                   ZWORK3D(:,:,:)=ZVARSAVE(:,:,:,J4,J5,J6)
1020                   ikdebzint=2
1021                   IF (INDEX(YTYPEOUT(1:4),'Z')/=0 .OR. INDEX(YTYPEOUT(1:4),'z')/=0) THEN
1022                     CALL ZINTER(ZWORK3D,XZZ,ZVARZCST,zlistevert,ikdebzint,XSPVAL)
1023                   ELSE IF (INDEX(YTYPEOUT(1:4),'P')/=0 .OR. INDEX(YTYPEOUT(1:4),'p')/=0) THEN
1024                     CALL PINTER(ZWORK3D,IGRID,XSPVAL,zlistevert,ZVARZCST,ZPABS)
1025                   ELSE IF (INDEX(YTYPEOUT(1:4),'H')/=0 .OR. INDEX(YTYPEOUT(1:4),'h')/=0) THEN
1026                     ZVARZCST(:,:,:)=ZWORK3D(:,:,:)
1027                   ELSE
1028                     print*,'** ERROR in vertical interpolations with ',YTYPEOUT
1029                   ENDIF
1030                   XVAR(:,:,:,J4,J5,J6)=ZVARZCST
1031                 END DO
1032               END DO
1033             END DO
1034             DEALLOCATE(ZVARSAVE,ZVARZCST,ZWORK3D)
1035             zmini=MINVAL(XVAR(:,:,:,:,:,:),MASK=XVAR(:,:,:,:,:,:)/=XSPVAL)
1036             zmaxi=MAXVAL(XVAR(:,:,:,:,:,:),MASK=XVAR(:,:,:,:,:,:)/=XSPVAL)
1037             print * ,' After vertical interpolation, min,max of the variable ',TRIM(CGROUP),'=', zmini,zmaxi
1038             ivarkdeb=1
1039             ivarkfin=inbvertz
1040             IF (ilocverbia >= 5 ) then
1041               print*,'ivarkdeb,ivarkfin= ',ivarkdeb,ivarkfin 
1042             ENDIF
1043           ENDIF
1044         ENDIF
1045         ! c. interpolation eventuelle sur l horizontale
1046         IF ( YOUTGRID(1:4) == 'LALO' ) THEN
1047           if (ilocverbia >= 0 ) then
1048             print *,'Translate to a regular lat lon grid '
1049           end if
1050           IF ( .NOT. ALLOCATED (ZNEWX) ) THEN
1051             IF ( IFLAGzcst == 1 ) THEN
1052               print*,'** no processing with ALT at the first group'
1053               GOTO 99
1054             ELSE
1055             ! c.1. creation de la grille réguliere en lat lon
1056               if (ilocverbia >= 2 ) then
1057                 print *,'grid creation, size of XLON: ',SIZE(XLON,1),SIZE(XLON,2) 
1058               end if
1059               ! calcul des coord X Y des points de la grille lat-lon reguliere
1060               ! determine le maximum d espacement en lat et lon sur le zoom
1061               ZDELTALON=max(XLON(ivarideb+1,ivarjdeb)-XLON(ivarideb,ivarjdeb)&
1062                            ,XLON(ivarifin,ivarjfin)-XLON(ivarifin-1,ivarjfin))
1063               ZDELTALAT=max(XLAT(ivarideb,ivarjdeb+1)-XLAT(ivarideb,ivarjdeb)&
1064                            ,XLAT(ivarifin,ivarjfin)-XLAT(ivarifin,ivarjfin-1))
1065               if (ZDELTALON == 0 .OR. ZDELTALAT == 0 ) THEN
1066                 print *,' error during ZDELTALON,ZDELTALAT computation=', ZDELTALON,ZDELTALAT
1067                 print *,'XLON(ivarideb+1,ivarjdeb)-XLON(ivarideb,ivarjdeb)'&
1068                         ,'XLON(ivarifin,ivarjfin)-XLON(ivarifin-1,ivarjfin)'&
1069                         ,'XLAT(ivarideb,ivarjdeb+1)-XLAT(ivarideb,ivarjdeb)'&
1070                         ,'XLAT(ivarifin,ivarjfin)-XLAT(ivarifin,ivarjfin-1)'
1071                 print *,XLON(ivarideb+1,ivarjdeb)-XLON(ivarideb,ivarjdeb)&
1072                        ,XLON(ivarifin,ivarjfin)-XLON(ivarifin-1,ivarjfin)&
1073                        ,XLAT(ivarideb,ivarjdeb+1)-XLAT(ivarideb,ivarjdeb)&
1074                        ,XLAT(ivarifin,ivarjfin)-XLAT(ivarifin,ivarjfin-1)           
1075                 print *, 'ivarideb+1,ivarjdeb,ivarifin-1,ivarjfin',ivarideb+1,ivarjdeb,ivarifin-1,ivarjfin
1076                 print *,'Verify the fields LAT LON of the FM file'
1077                 ALLOCATE(ZX(SIZE(XLAT,1),SIZE(XLAT,2)),ZY(SIZE(XLAT,1),SIZE(XLAT,2)))
1078                 ZX(1:SIZE(XZZ,1),1) = XXX(1:SIZE(XZZ,1),IGRID)
1079                 ZX(:,2:SIZE(XZZ,2)) = SPREAD(ZX(:,1),2,SIZE(XZZ,2)-1)
1080                 ZY(1,1:SIZE(XZZ,2)) = XXY(1:SIZE(XZZ,2),IGRID)
1081                 ZY(2:SIZE(XZZ,1),:) = SPREAD(ZY(1,:),1,SIZE(XZZ,1)-1)
1082                 !CALL SM_LATLON(XXHAT,XYHAT,XLATORI,XLONORI, &
1083                 !! XXHAT,XYHAT supprimes en masdev4_7
1084                 CALL SM_LATLON(XLATORI,XLONORI,ZX,ZY,XLAT,XLON)
1085                 ZDELTALON=max(XLON(ivarideb+1,ivarjdeb)-XLON(ivarideb,ivarjdeb)&
1086                            ,XLON(ivarifin,ivarjfin)-XLON(ivarifin-1,ivarjfin))
1087                 ZDELTALAT=max(XLAT(ivarideb,ivarjdeb+1)-XLAT(ivarideb,ivarjdeb)&
1088                            ,XLAT(ivarifin,ivarjfin)-XLAT(ivarifin,ivarjfin-1))
1089                 print *,' After Model Grid computation: ZDELTALON,ZDELTALAT=', ZDELTALON,ZDELTALAT
1090               endif
1091               IDIM1=(maxval(XLON)-minval(XLON))/ZDELTALON
1092               IDIM2=(maxval(XLAT)-minval(XLAT))/ZDELTALAT
1093               ALLOCATE (ZNEWLAT(IDIM1,IDIM2),ZNEWLON(IDIM1,IDIM2) )
1094               if (ilocverbia >= 1 ) then
1095                 print*,' ZDELTALON,ZDELTALAT= ',ZDELTALON,ZDELTALAT
1096               endif
1097               if (ilocverbia >= 2 ) then
1098                 print*,' IDIM1,IDIM2= ',IDIM1,IDIM2
1099               endif
1100               ! depart de la nouvelle grille : coin Sud Ouest
1101               DO JI=1,IDIM1
1102                 ZNEWLON(JI,:)=minval(XLON) + (JI-1) *ZDELTALON
1103               ENDDO
1104               DO JJ=1,IDIM2
1105                 ZNEWLAT(:,JJ)=minval(XLAT) + (JJ-1) *ZDELTALAT
1106               ENDDO
1107               if (ilocverbia >= 4 ) then
1108                 print*, 'new lat lon grid=',ZNEWLAT(1,:)
1109                 print*, ZNEWLON(:,1)
1110               endif
1111               ALLOCATE (ZNEWX(IDIM1,IDIM2))
1112               ALLOCATE (ZNEWY(IDIM1,IDIM2))
1113               CALL SM_XYHAT(XLATORI,XLONORI,ZNEWLAT,ZNEWLON,ZNEWX,ZNEWY)
1114               if (ilocverbia >= 4 ) then
1115                 ! XXX= XXHAT et XXY=XYHAT pour les 7 grilles
1116                 print*,' After SM_XYHAT old limits X: ', &
1117                        XXX(1,IGRID),XXX(SIZE(XVAR,1),IGRID)
1118                 print*,'                new limits X: ', &
1119                        ZNEWX(1,1),ZNEWX(IDIM1,IDIM2)
1120                 print*,'                old limits Y: ', &
1121                        XXY(1,IGRID),XXY(SIZE(XVAR,2),IGRID)
1122                 print*,'                new limits Y: ', &
1123                        ZNEWY(1,1),ZNEWY(IDIM1,IDIM2)
1124               endif
1125               if (ilocverbia >= 5 ) then
1126                 DO JI= 1,SIZE(XVAR,1) 
1127                   print*,'XXHAT ZNEWX',XXX(JI,IGRID),ZNEWX(JI,1),ZNEWX(JI,IDIM2)
1128                 ENDDO
1129                 DO JJ= 1,SIZE(XVAR,2) 
1130                   print*,'XYHAT ZNEWY',XXY(JJ,IGRID),ZNEWY(1,JJ),ZNEWX(IDIM1,JJ)
1131                 ENDDO
1132               endif
1133               ! calcul de la section de tableau correspondant au zoom
1134               I1=(maxval(XLON(ivarideb:ivarifin,ivarjdeb:ivarjfin)) &
1135                  -minval(XLON(ivarideb:ivarifin,ivarjdeb:ivarjfin)) )/ZDELTALON
1136               I2=(maxval(XLAT(ivarideb:ivarifin,ivarjdeb:ivarjfin)) &
1137                  -minval(XLAT(ivarideb:ivarifin,ivarjdeb:ivarjfin)) )/ZDELTALAT
1138               IZOOMIDEB=MAX(MIN(COUNT(ZNEWLON(:,1)<XLON(ivarideb,ivarjdeb)),IDIM1),1)
1139               IZOOMJDEB=MAX(MIN(COUNT(ZNEWLAT(1,:)<XLAT(ivarideb,ivarjdeb)),IDIM2),1)
1140               !IZOOMIFIN=MIN(IZOOMIDEB+I1,IDIM1)
1141               !IZOOMJFIN=MIN(IZOOMJDEB+I2,IDIM2)
1142               IZOOMIFIN=MAX(MIN(COUNT(ZNEWLON(:,1)<XLON(ivarifin,ivarjfin)),IDIM1),1)
1143               IZOOMJFIN=MAX(MIN(COUNT(ZNEWLAT(1,:)<XLAT(ivarifin,ivarjfin)),IDIM2),1)
1144               if (ilocverbia >= 2 ) then
1145                 print*,' ZOOM along i in the LON-LAT grid: ', &
1146                        IZOOMIDEB,IZOOMIFIN,I1
1147                 print*,'            j                    : ', &
1148                        IZOOMJDEB,IZOOMJFIN,I2
1149               endif
1150             ENDIF
1151           ENDIF ! fin grille ZNEWX deja allouee
1152           ! c.2. interpolation sur la nouvelle grille
1153           IF( IFLAGzcst/= 1 .AND. (NREADIH-NREADIL)>0 .AND. (NREADJH-NREADJL)>0 )THEN
1154             ! interpolation vers la nouvelle grille réguliere en lat lon 
1155             !sauf la grille verticale definie en niveaux Z et champs 1D
1156             if (ilocverbia >= 1 ) then
1157               print*,' interpolation for the variable  ',trim(CGROUP)
1158             end if
1159             allocate(ZWORK3D(SIZE(XVAR,1),SIZE(XVAR,2),SIZE(XVAR,3)))
1160             allocate(ZWORK3D2(IDIM1,IDIM2,SIZE(XVAR,3)))
1161             ! stockage des champs interpoles dans la nouvelle grille
1162             if (allocated (ZVARSAVE)) DEALLOCATE(ZVARSAVE)
1163             allocate(ZVARSAVE(IDIM1,IDIM2,SIZE(XVAR,3),&
1164                      SIZE(XVAR,4),SIZE(XVAR,5),SIZE(XVAR,6)))
1165             ! boucle sur les dimensions 4 5 6
1166             DO J6=ivarprocinf,ivarprocsup
1167               DO J5=ivartrajinf,ivartrajsup
1168                 DO J4=ivartinf,ivartsup
1169                   ZWORK3D(:,:,:)=XVAR(:,:,:,J4,J5,J6)
1170                   if (ilocverbia >= 2 ) then
1171                     print *,'before HOR_INTERP_4PTS J4,J5,J6=', J4,J5,J6
1172                   end if
1173                   CALL HOR_INTERP_4PTS(XXX(:,IGRID),XXY(:,IGRID),ZWORK3D, &
1174                                        ZNEWX,ZNEWY,ZWORK3D2)
1175                   ZVARSAVE(:,:,:,J4,J5,J6)=ZWORK3D2(:,:,:)
1176                 END DO
1177               END DO
1178             END DO
1179             ! resultat dans XVAR passe en module
1180             DEALLOCATE (XVAR)
1181             ALLOCATE(XVAR(IDIM1,IDIM2,SIZE(ZVARSAVE,3),&
1182                     SIZE(ZVARSAVE,4),SIZE(ZVARSAVE,5),SIZE(ZVARSAVE,6)))
1183             XVAR=XSPVAL
1184             XVAR(:,:,:,ivartinf:ivartsup,ivartrajinf:ivartrajsup,ivarprocinf:ivarprocsup)= &
1185             ZVARSAVE(:,:,:,ivartinf:ivartsup,ivartrajinf:ivartrajsup,ivarprocinf:ivarprocsup)
1186             DEALLOCATE (ZVARSAVE)
1187             zmini=MINVAL(XVAR(:,:,:,:,:,:),MASK=XVAR(:,:,:,:,:,:)/=XSPVAL)
1188             zmaxi=MAXVAL(XVAR(:,:,:,:,:,:),MASK=XVAR(:,:,:,:,:,:)/=XSPVAL)
1189             print * ,' After horizontal interpolation, min,max of the variable ',TRIM(CGROUP),'=', zmini,zmaxi
1190             if (ilocverbia >= 2 ) then
1191               print*, 'After HOR_INTERP_4PTS all the dim 4,5,6'
1192             endif
1193             IF (allocated(ZWORK3D)) DEALLOCATE(ZWORK3D)
1194             IF (allocated(ZWORK3D2)) DEALLOCATE(ZWORK3D2)
1195           ENDIF
1196         ENDIF
1197         ! d. ecriture des donnees au format cdl ou llz/llp
1198         IF ( YOUTGRID(1:4) == 'LALO' ) THEN
1199           IF ( IFLAGzcst /= 1 ) THEN
1200             ivarideb=IZOOMIDEB
1201             ivarifin=IZOOMIFIN
1202             ivarjdeb=IZOOMJDEB
1203             ivarjfin=IZOOMJFIN
1204           ENDIF
1205           SELECT CASE(YTYPEOUT(1:4))
1206           CASE('LLZV','llzv','LLPV','llpv')
1207             IF (allocated(ZWORK3D)) DEALLOCATE(ZWORK3D)
1208             ALLOCATE(ZWORK3D(size(XVAR,1),size(XVAR,2),size(XVAR,3)))
1209             IF (SIZE(XVAR,3)==inbvertz) THEN
1210               ZWORK3D(1,1,:)=zlistevert
1211             ELSE
1212               ZWORK3D(1,1,:)=XSPVAL
1213             ENDIF
1214             CALL WRITELLHV(ivarideb,ivarifin,ivarjdeb,ivarjfin, &
1215                        ivarkdeb,ivarkfin,ivartinf,ivartsup, &
1216                        ivartrajinf,ivartrajsup,ivarprocinf,ivarprocsup,&
1217                        CGROUP,YFILEIN,YFLAGWRITE,YTYPEOUT,&
1218                        ilocverbia,iret,PLON=ZNEWLON,PLAT=ZNEWLAT,&
1219                        PALT=ZWORK3D)       
1220             if (ilocverbia > 0 ) then
1221               print*,'WRITELLHV LALO return= ', YTYPEOUT,'= ',iret
1222             end if
1223             DEALLOCATE(ZWORK3D)
1224         !
1225           CASE('KCDL','ZCDL','PCDL')
1226            YGROUP=ADJUSTL(ADJUSTR(CGROUP)//ADJUSTL(YK))
1227            CALL WRITECDL(ivarideb,ivarifin,ivarjdeb,ivarjfin, &
1228                         ivarkdeb,ivarkfin,ivartinf,ivartsup, &
1229                         ivartrajinf,ivartrajsup,ivarprocinf,ivarprocsup, &
1230                         YGROUP,YFILEIN,YFLAGWRITE,YOUTGRID,YSUFFIX_file, &
1231                         ilocverbia,iret,PGRIDX=ZNEWLON(:,1),PGRIDY=ZNEWLAT(1,:))
1232            IF (ilocverbia >= 1 ) print *,' counter of added fields=',inetadd
1233            if ( inetadd == 0) then
1234              print *,' The program adds the ALT 3Dfield to the netcdf file'
1235              YGROUP_OLD=CGROUP(1:13)
1236              CGROUP='ALT'
1237              inetadd=inetadd+1
1238              YFLAGWRITE='OLD'
1239              ino_init_zoom=1
1240              GO TO 77
1241            endif
1242            if ( inetadd == 1 .AND. YOUTGRID(1:4) == 'CONF' )THEN
1243              print *,' The program adds the LAT 3Dfield to the netcdf file'      
1244              CGROUP='LAT'
1245              inetadd=inetadd+1
1246              ino_init_zoom=1
1247              GO TO 77
1248            endif
1249            if ( inetadd == 2 .AND. YOUTGRID(1:4) == 'CONF' )THEN
1250              print *,' The program adds the LON 3Dfield to the netcdf file'      
1251              CGROUP='LON'
1252              inetadd=inetadd+1
1253              ino_init_zoom=1
1254              GO TO 77
1255            endif
1256            
1257           END SELECT
1258         ELSE ! pas d interpolation horizontale
1259           SELECT CASE(YTYPEOUT(1:4))
1260           CASE('LLZV','llzv','LLPV','llpv')
1261             IF (SIZE(XVAR,3)==inbvertz) THEN  ! champ 3D
1262               IF (allocated(ZWORK3D)) DEALLOCATE(ZWORK3D)
1263               ALLOCATE(ZWORK3D(size(XVAR,1),size(XVAR,2),size(XVAR,3)))
1264               ZWORK3D(1,1,:)=zlistevert
1265               CALL WRITELLHV(ivarideb,ivarifin,ivarjdeb,ivarjfin, &
1266                        ivarkdeb,ivarkfin,ivartinf,ivartsup, &
1267                        ivartrajinf,ivartrajsup,ivarprocinf,ivarprocsup,&
1268                        CGROUP,YFILEIN,YFLAGWRITE,YTYPEOUT,&
1269                        ilocverbia,iret,&
1270                        PALT=ZWORK3D)       
1271             ELSE                              ! champ 2D
1272               IF((YTYPEOUT(3:3)=='z').OR.(YTYPEOUT(3:3)=='p')) YTYPEOUT3='h'
1273               IF((YTYPEOUT(3:3)=='Z').OR.(YTYPEOUT(3:3)=='P')) YTYPEOUT3='H'
1274               CALL WRITELLHV(ivarideb,ivarifin,ivarjdeb,ivarjfin, &
1275                        ivarkdeb,ivarkfin,ivartinf,ivartsup, &
1276                        ivartrajinf,ivartrajsup,ivarprocinf,ivarprocsup, &
1277                        CGROUP,YFILEIN,YFLAGWRITE, &
1278                        YTYPEOUT(1:2)//YTYPEOUT3//YTYPEOUT(4:4), &
1279                        ilocverbia,iret)
1280             ENDIF
1281             if (ilocverbia > 0 ) then
1282               print*,' WRITELLHV for ', YTYPEOUT,', return value= ',iret
1283             end if
1284             IF (allocated(ZWORK3D)) DEALLOCATE(ZWORK3D)
1285         !
1286           CASE('KCDL','ZCDL','PCDL')
1287             YGROUP=ADJUSTL(ADJUSTR(CGROUP)//ADJUSTL(YK))
1288             CALL WRITECDL(ivarideb,ivarifin,ivarjdeb,ivarjfin, &
1289                         ivarkdeb,ivarkfin,ivartinf,ivartsup, &
1290                         ivartrajinf,ivartrajsup,ivarprocinf,ivarprocsup, &
1291                         YGROUP,YFILEIN,YFLAGWRITE,YOUTGRID,YSUFFIX_file, &
1292                         ilocverbia,iret,PGRIDX=XXX(:,IGRID),PGRIDY=XXY(:,IGRID))
1293             IF (ilocverbia >= 1 ) print *,' counter of added fields=',inetadd
1294            if ( inetadd == 0) then
1295              if (ivarkdeb == ivarkfin .AND. ivarkdeb == 1 ) THEN
1296                print *, 'No ALT field for only one vertical position'
1297              else
1298              print *,' The program adds the ALT 3Dfield to the netcdf file'
1299              YGROUP_OLD=CGROUP(1:13)
1300              CGROUP='ALT'
1301              inetadd=inetadd+1
1302              YFLAGWRITE='OLD'
1303              ino_init_zoom=1
1304              GO TO 77
1305              endif
1306            endif
1307            if ( inetadd == 1 .AND. YOUTGRID(1:4) == 'CONF' )THEN
1308              if (ivarideb /= ivarifin .AND. ivarjdeb /= ivarjfin ) THEN
1309
1310              print *,' The program adds the LAT 3Dfield to the netcdf file'      
1311              CGROUP='LAT'
1312              inetadd=inetadd+1
1313              ino_init_zoom=1
1314              GO TO 77
1315              else
1316                print *, ' No LAT field for only one location', ivarideb,ivarifin,ivarjdeb,ivarjfin
1317              endif
1318            endif
1319            if ( inetadd == 2 .AND. YOUTGRID(1:4) == 'CONF' )THEN
1320              if (ivarideb /= ivarifin .AND. ivarjdeb /= ivarjfin ) THEN      
1321              print *,' The program adds the LON 3Dfield to the netcdf file'      
1322              CGROUP='LON'
1323              inetadd=inetadd+1
1324              ino_init_zoom=1
1325              GO TO 77
1326              else
1327                print *, ' No LON field for only one location', ivarideb,ivarifin,ivarjdeb,ivarjfin
1328              endif
1329            endif
1330           END SELECT              
1331         ENDIF
1332         ! retour a XZZ pour NGRID a 4 (cf readvar)
1333         CALL COMPCOORD_FORDIACHRO(4)
1334      END SELECT
1335      ! indiquera aux routines d ecriture que le fichier courant est deja ouvert
1336      YFLAGWRITE='OLD'
1337   ! 
1338   ELSE   ! iret /=0
1339     print *, ' READVAR return= ',iret
1340   ENDIF  
1341 END DO ! boucle champ a traiter
1342 !
1343 !
1344 !---------------------------------------------------------------------------
1345 !
1346 !*       4.    CLOSURE OF OUTPUT FILE
1347 !              ----------------------
1348 !
1349 !pour clore le traitement meme si la liste des champs est non terminee par END
1350 88 CONTINUE
1351 !
1352 IF (ALLOCATED(ZNEWX))   DEALLOCATE(ZNEWX,ZNEWY)
1353 IF (ALLOCATED(ZNEWLAT)) DEALLOCATE(ZNEWLAT,ZNEWLON)
1354 IF (ALLOCATED(ZWORK2D)) DEALLOCATE(ZWORK2D,ZWORK2D2)
1355 !
1356 PRINT*, 'END ->  Close the output file'
1357 YFLAGWRITE='CLO'
1358 SELECT CASE(YTYPEOUT(1:4))
1359   CASE('DIAC')
1360     CALL WRITEVAR(ivarideb,ivarifin,ivarjdeb,ivarjfin,ivarkdeb,ivarkfin,&
1361                   ivartinf,ivartsup,ivartrajinf,ivartrajsup,ivarprocinf,ivarprocsup,  &
1362                  CGROUP,YFILEIN,YFLAGWRITE,'2  ',ilocverbia,iret)
1363   CASE('LLHV','llhv','LLZV','llzv','LLPV','llpv')             
1364     CALL WRITELLHV(ivarideb,ivarifin,ivarjdeb,ivarjfin,ivarkdeb,ivarkfin,&
1365                  ivartinf,ivartsup,ivartrajinf,ivartrajsup,ivarprocinf,ivarprocsup,&
1366                  CGROUP,YFILEIN,YFLAGWRITE,YTYPEOUT,ilocverbia,iret)      
1367   CASE('KCDL','ZCDL','PCDL')
1368     CALL WRITECDL(ivarideb,ivarifin,ivarjdeb,ivarjfin,ivarkdeb,ivarkfin,&
1369                   ivartinf,ivartsup,ivartrajinf,ivartrajsup,ivarprocinf,ivarprocsup,&
1370                   CGROUP,YFILEIN,YFLAGWRITE,YOUTGRID,YSUFFIX_file,      &
1371                   ilocverbia,iret,PGRIDX=XXX(:,IGRID),PGRIDY=XXY(:,IGRID))
1372   CASE DEFAULT
1373     PRINT*, 'Closure of output type ',YTYPEOUT ,' not coded'
1374 END SELECT
1375 !
1376 !-------------------------------------------------------------------------------
1377 !
1378 !*       5.    END
1379 !              ---
1380 !
1381 99 CONTINUE
1382 PRINT*, 'Delete the links if necessary'
1383 YDUMMYFILE=''
1384 CALL CREATLINK(' ',YDUMMYFILE,'CLEAN',ILOCVERBIA)
1385 PRINT*, 'The file ',TRIM(YLUDIR),' stores all the input directives '
1386 PRINT*, ' you must give a new name to use it again'
1387 CLOSE(ILUDIR)
1388 !
1389 !-------------------------------------------------------------------------------
1390 !
1391 END PROGRAM EXTRACTDIA
1392 !