ab29f4427e96f270667bbb8fca9f929f9732dce3
[MNH-git_open_source-lfs.git] / tools / lfi2cdf / src / mode_util.f90
1 MODULE mode_util
2   USE MODE_FIELDTYPE
3   USE mode_dimlist
4   USE MODD_PARAM
5   USE netcdf
6
7   IMPLICIT NONE 
8
9   TYPE workfield
10      CHARACTER(LEN=FM_FIELD_SIZE)            :: name   ! nom du champ
11      INTEGER                                 :: TYPE   ! type (entier ou reel)    
12      CHARACTER(LEN=:), POINTER               :: comment
13      TYPE(dimCDF),                   POINTER :: dim
14      INTEGER                                 :: id
15      INTEGER                                 :: grid
16   END TYPE workfield
17
18 #ifndef LOWMEM
19   TYPE lfidata
20      INTEGER(KIND=8), DIMENSION(:), POINTER :: iwtab
21   END TYPE lfidata
22   TYPE(lfidata), DIMENSION(:), ALLOCATABLE :: lfiart
23 #endif
24
25   LOGICAL(KIND=LFI_INT), PARAMETER :: ltrue  = .TRUE.
26   LOGICAL(KIND=LFI_INT), PARAMETER :: lfalse = .FALSE.
27
28 CONTAINS 
29   FUNCTION str_replace(hstr, hold, hnew)
30     CHARACTER(LEN=*) :: hstr, hold, hnew
31     CHARACTER(LEN=LEN_TRIM(hstr)+MAX(0,LEN(hnew)-LEN(hold))) :: str_replace
32     
33     INTEGER :: pos
34     
35     pos = INDEX(hstr,hold)
36     IF (pos /= 0) THEN
37        str_replace = hstr(1:pos-1)//hnew//hstr(pos+LEN(hold):)
38     ELSE 
39        str_replace = hstr 
40     END IF
41
42   END FUNCTION str_replace
43
44   SUBROUTINE FMREADLFIN1(klu,hrecfm,kval,kresp)
45   INTEGER(KIND=LFI_INT), INTENT(IN) :: klu ! logical fortran unit au lfi file
46   CHARACTER(LEN=*),INTENT(IN)       :: hrecfm ! article name to be read
47   INTEGER, INTENT(OUT)        :: kval ! integer value for hrecfm article
48   INTEGER(KIND=LFI_INT), INTENT(OUT):: kresp! return code null if OK
49   !
50   INTEGER(KIND=8),DIMENSION(:),ALLOCATABLE::iwork
51   INTEGER :: icomlen
52   INTEGER(KIND=LFI_INT) :: iresp,ilenga,iposex
53   !
54   CALL LFINFO(iresp,klu,hrecfm,ilenga,iposex)
55   IF (iresp /=0 .OR. ilenga == 0) THEN
56     kresp = -1
57     kval = 0
58   ELSE
59     ALLOCATE(IWORK(ilenga))
60     CALL LFILEC(iresp,klu,hrecfm,iwork,ilenga)
61     icomlen = iwork(2)
62     kval = iwork(3+icomlen)
63     kresp = iresp
64     DEALLOCATE(IWORK)
65   END IF
66   END SUBROUTINE FMREADLFIN1
67
68   SUBROUTINE parse_lfi(klu, hvarlist, knaf, tpreclist, kbuflen, current_level)
69     INTEGER, INTENT(IN)                    :: klu
70     INTEGER, INTENT(INOUT)                 :: knaf
71     CHARACTER(LEN=*), intent(IN)           :: hvarlist
72     TYPE(workfield), DIMENSION(:), POINTER :: tpreclist    
73     INTEGER, INTENT(OUT)                   :: kbuflen
74     INTEGER, INTENT(IN), OPTIONAL          :: current_level
75
76     INTEGER                                  :: ji,jj
77     INTEGER                                  :: ndb, nde
78     INTEGER                                  :: inaf
79     LOGICAL                                  :: ladvan
80     INTEGER                                  :: ich
81     INTEGER                                  :: fsize,sizemax
82     CHARACTER(LEN=FM_FIELD_SIZE)             :: yrecfm
83     CHARACTER(LEN=4)                         :: suffix
84 #ifdef LOWMEM
85     INTEGER(KIND=8),DIMENSION(:),ALLOCATABLE :: iwork
86 #endif
87     INTEGER(KIND=LFI_INT)                    :: iresp,ilu,ileng,ipos
88     !JUAN CYCCL3
89     INTEGER                        :: JPHEXT
90
91     ilu = klu
92
93     CALL FMREADLFIN1(klu,'JPHEXT',JPHEXT,iresp)
94     IF (iresp /= 0) JPHEXT=1
95     ! First check if IMAX,JMAX,KMAX exist in LFI file
96     ! to handle 3D, 2D variables -> update IDIMX,IDIMY,IDIMZ
97     CALL FMREADLFIN1(klu,'IMAX',IDIMX,iresp)
98     IF (iresp == 0) IDIMX = IDIMX+2*JPHEXT  ! IMAX + 2*JPHEXT
99     !
100     CALL FMREADLFIN1(klu,'JMAX',IDIMY,iresp)
101     IF (iresp == 0) IDIMY = IDIMY+2*JPHEXT  ! JMAX + 2*JPHEXT
102     !
103     CALL FMREADLFIN1(ilu,'KMAX',IDIMZ,iresp)
104     IF (iresp == 0) IDIMZ = IDIMZ+2  ! KMAX + 2*JPVEXT
105     GUSEDIM = (IDIMX*IDIMY > 0)
106     IF (GUSEDIM) THEN
107       PRINT *,'MESONH 3D, 2D articles DIMENSIONS used :'
108       PRINT *,'DIMX =',IDIMX
109       PRINT *,'DIMY =',IDIMY
110       PRINT *,'DIMZ =',IDIMZ ! IDIMZ may be equal to 0 (PGD files)
111     ELSE
112       PRINT *,'BEWARE : ALL MesoNH arrays are handled as 1D arrays !'
113     END IF
114
115 #ifndef LOWMEM
116     ALLOCATE(lfiart(knaf))
117 #endif
118     ALLOCATE(tpreclist(knaf))
119     sizemax = 0
120
121     IF (present(current_level)) THEN
122       write(suffix,'(I4.4)') current_level
123     ElSE
124       suffix=''
125     END IF
126
127     ! Phase 1 : build articles list to convert.
128     !
129     !    Pour l'instant tous les articles du fichier LFI sont
130     !    convertis. On peut modifier cette phase pour prendre en
131     !    compte un sous-ensemble d'article (liste definie par
132     !    l'utilisateur par exemple)  
133     !
134     IF (LEN_TRIM(hvarlist) > 0) THEN
135        ! A variable list is provided with -v var1,...
136        ndb  = 1
137        inaf = 0
138        DO ji=1,knaf
139           nde = INDEX(TRIM(hvarlist(ndb:)),',')
140           yrecfm = hvarlist(ndb:ndb+nde-2)
141           ndb = nde+ndb
142
143           CALL LFINFO(iresp,ilu,trim(yrecfm)//trim(suffix),ileng,ipos)
144           
145           IF (iresp /= 0 .OR. ileng == 0) THEN
146              PRINT *,'Article ',TRIM(yrecfm), ' not found!'
147           ELSE
148              inaf = inaf+1
149              ! PRINT *,'Article ',ji,' : ',TRIM(yrecfm),', longueur = ',ileng
150              tpreclist(inaf)%name = yrecfm
151              IF (ileng > sizemax) sizemax = ileng        
152 #ifndef LOWMEM       
153              ALLOCATE(lfiart(inaf)%iwtab(ileng))
154 #endif
155           end IF
156        END DO
157     ELSE
158        ! Entire file is converted
159        CALL LFIPOS(iresp,ilu)
160        ladvan = .TRUE.
161        
162        DO ji=1,knaf
163           CALL LFICAS(iresp,ilu,yrecfm,ileng,ipos,ladvan)
164           ! PRINT *,'Article ',ji,' : ',TRIM(yrecfm),', longueur = ',ileng
165           tpreclist(ji)%name = yrecfm
166           IF (ileng > sizemax) sizemax = ileng        
167 #ifndef LOWMEM       
168           ALLOCATE(lfiart(ji)%iwtab(ileng))
169 #endif
170        END DO
171        inaf = knaf
172     END IF
173
174     kbuflen = sizemax
175 #ifdef LOWMEM
176     WRITE(*,'("Taille maximale du buffer :",f10.3," Mo")') sizemax*8./1048576.
177     ALLOCATE(iwork(sizemax))
178 #endif
179     
180     ! Phase 2 : Extract comments and dimensions for valid articles.
181     !           Infos are put in tpreclist.
182     CALL init_dimCDF()
183     DO ji=1,inaf
184        yrecfm = trim(tpreclist(ji)%name)//trim(suffix)
185        CALL LFINFO(iresp,ilu,yrecfm,ileng,ipos)
186 #ifdef LOWMEM
187        CALL LFILEC(iresp,ilu,yrecfm,iwork,ileng)
188        tpreclist(ji)%TYPE = get_ftype(yrecfm)               
189        tpreclist(ji)%grid = iwork(1)
190
191        ALLOCATE(character(len=iwork(2)) :: tpreclist(ji)%comment)
192        DO jj=1,iwork(2)
193           ich = iwork(2+jj)
194           tpreclist(ji)%comment(jj:jj) = CHAR(ich)
195        END DO
196        fsize = ileng-(2+iwork(2))
197 #else
198        CALL LFILEC(iresp,ilu,yrecfm,lfiart(ji)%iwtab,ileng)
199        tpreclist(ji)%TYPE = get_ftype(yrecfm)               
200        tpreclist(ji)%grid = lfiart(ji)%iwtab(1)
201
202        ALLOCATE(character(len=lfiart(ji)%iwtab(2)) :: tpreclist(ji)%comment)
203        DO jj=1,lfiart(ji)%iwtab(2)
204           ich = lfiart(ji)%iwtab(2+jj)
205           tpreclist(ji)%comment(jj:jj) = CHAR(ich)
206        END DO
207        fsize = ileng-(2+lfiart(ji)%iwtab(2))
208 #endif
209        tpreclist(ji)%dim=>get_dimCDF(fsize)
210     END DO
211   
212     PRINT *,'Nombre de dimensions = ', size_dimCDF()
213 #ifdef LOWMEM
214     DEALLOCATE(iwork)
215 #endif
216     knaf = inaf
217   END SUBROUTINE parse_lfi
218   
219   SUBROUTINE read_data_lfi(klu, hvarlist, nbvar, tpreclist, kbuflen, current_level)
220     INTEGER, INTENT(IN)                    :: klu
221     INTEGER, INTENT(INOUT)                 :: nbvar
222     CHARACTER(LEN=*), intent(IN)           :: hvarlist
223     TYPE(workfield), DIMENSION(:), POINTER :: tpreclist
224     INTEGER, INTENT(IN)                    :: kbuflen
225     INTEGER, INTENT(IN), OPTIONAL          :: current_level
226
227     INTEGER                                  :: ji,jj
228     INTEGER                                  :: ndb, nde
229     LOGICAL                                  :: ladvan
230     INTEGER                                  :: ich
231     INTEGER                                  :: fsize,sizemax
232     CHARACTER(LEN=FM_FIELD_SIZE)             :: yrecfm
233     CHARACTER(LEN=4)                         :: suffix
234 #ifdef LOWMEM
235     INTEGER(KIND=8),DIMENSION(:),ALLOCATABLE :: iwork
236 #endif
237     INTEGER(KIND=LFI_INT)                    :: iresp,ilu,ileng,ipos
238
239     ilu = klu
240
241     IF (present(current_level)) THEN
242       write(suffix,'(I4.4)') current_level
243     ElSE
244       suffix=''
245     END IF
246
247 #ifdef LOWMEM
248     ALLOCATE(iwork(kbuflen))
249 #endif
250
251     DO ji=1,nbvar
252        yrecfm = trim(tpreclist(ji)%name)//trim(suffix)
253        CALL LFINFO(iresp,ilu,yrecfm,ileng,ipos)
254 #ifdef LOWMEM
255        CALL LFILEC(iresp,ilu,yrecfm,iwork,ileng)
256        tpreclist(ji)%grid = iwork(1)
257 #else
258        CALL LFILEC(iresp,ilu,yrecfm,lfiart(ji)%iwtab,ileng)
259        tpreclist(ji)%grid = lfiart(ji)%iwtab(1)
260 #endif
261     END DO
262
263 #ifdef LOWMEM
264     DEALLOCATE(iwork)
265 #endif
266   END SUBROUTINE read_data_lfi
267
268   SUBROUTINE HANDLE_ERR(status,line)
269     INTEGER :: status,line
270
271     IF (status /= NF90_NOERR) THEN
272        PRINT *, 'line ',line,': ',NF90_STRERROR(status)
273        STOP
274     END IF
275   END SUBROUTINE HANDLE_ERR
276
277   SUBROUTINE def_ncdf(tpreclist,knaf,oreduceprecision,kcdf_id,omerge,ocompress,compress_level)
278     TYPE(workfield),DIMENSION(:),INTENT(INOUT) :: tpreclist
279     INTEGER,                     INTENT(IN) :: knaf
280     LOGICAL,                     INTENT(IN) :: oreduceprecision
281     INTEGER,                     INTENT(OUT):: kcdf_id
282     LOGICAL,                     INTENT(IN) :: omerge
283     LOGICAL,                     INTENT(IN) :: ocompress
284     INTEGER,                     INTENT(IN) :: compress_level
285
286     INTEGER :: status
287     INTEGER :: ji
288     TYPE(dimCDF), POINTER :: tzdim
289     INTEGER               :: invdims
290     INTEGER               :: type_float
291     INTEGER, DIMENSION(10) :: ivdims
292     CHARACTER(LEN=20)     :: ycdfvar
293
294
295     IF (oreduceprecision) THEN
296       type_float = NF90_REAL
297     ELSE
298       type_float = NF90_DOUBLE
299     END IF
300
301       ! global attributes
302       status = NF90_PUT_ATT(kcdf_id,NF90_GLOBAL,'Title',VERSION_ID)
303       IF (status /= NF90_NOERR) CALL HANDLE_ERR(status,__LINE__)
304
305     ! define DIMENSIONS
306     tzdim=>first_DimCDF()
307     DO WHILE(ASSOCIATED(tzdim))
308       IF (tzdim%create) THEN
309         status = NF90_DEF_DIM(kcdf_id,tzdim%name,tzdim%len,tzdim%id)
310         IF (status /= NF90_NOERR) CALL HANDLE_ERR(status,__LINE__)
311       END IF
312       tzdim=>tzdim%next
313     END DO
314
315     PRINT *,'------------- NetCDF DEFINITION ---------------'
316
317     ! define VARIABLES and ATTRIBUTES
318     DO ji=1,knaf
319       
320        IF (ASSOCIATED(tpreclist(ji)%dim)) THEN
321          IF (tpreclist(ji)%dim%create) THEN
322            invdims   = 1
323            ivdims(1) = tpreclist(ji)%dim%id
324          ELSE
325            invdims = tpreclist(ji)%dim%ndims
326            IF(omerge) invdims=invdims+1 !when merging variables from LFI splitted files
327            SELECT CASE(invdims)
328            CASE(2)
329               ivdims(1)=ptdimx%id
330               ivdims(2)=ptdimy%id
331            CASE(3)
332               ivdims(1)=ptdimx%id
333               ivdims(2)=ptdimy%id
334               ivdims(3)=ptdimz%id
335            CASE(12)
336               ivdims(1)=ptdimx%id
337               ivdims(2)=ptdimz%id
338               invdims = 2 ! on retablit la bonne valeur du nbre de dimension
339            CASE default
340              PRINT *,'Fatal error in NetCDF dimension definition'
341              STOP
342            END SELECT
343          END IF
344        ELSE
345          ! scalar variables
346           invdims   = 0
347           ivdims(1) = 0 ! ignore dans ce cas
348        END IF
349        
350        ! Variables definition
351
352        !! NetCDF n'aime pas les '%' dans le nom des variables
353        !! "%" remplaces par '__' 
354        ycdfvar = str_replace(tpreclist(ji)%name,'%','__')
355        !! ni les '.' remplaces par '--'
356        ycdfvar = str_replace(ycdfvar,'.','--')
357
358        SELECT CASE(tpreclist(ji)%TYPE)
359        CASE (TEXT)
360 !          PRINT *,'TEXT : ',tpreclist(ji)%name
361           status = NF90_DEF_VAR(kcdf_id,ycdfvar,NF90_CHAR,&
362                    ivdims(:invdims),tpreclist(ji)%id)
363           IF (status /= NF90_NOERR) CALL HANDLE_ERR(status,__LINE__)
364
365        CASE (INT,BOOL)
366 !          PRINT *,'INT,BOOL : ',tpreclist(ji)%name
367           status = NF90_DEF_VAR(kcdf_id,ycdfvar,NF90_INT,&
368                    ivdims(:invdims),tpreclist(ji)%id)
369           IF (status /= NF90_NOERR) CALL HANDLE_ERR(status,__LINE__)
370
371        CASE(FLOAT)
372 !          PRINT *,'FLOAT : ',tpreclist(ji)%name
373           status = NF90_DEF_VAR(kcdf_id,ycdfvar,type_float,&
374                    ivdims(:invdims),tpreclist(ji)%id)
375           IF (status /= NF90_NOERR) CALL HANDLE_ERR(status,__LINE__)
376
377           
378        CASE default
379           PRINT *,'ATTENTION : ',TRIM(tpreclist(ji)%name),' est de&
380                & TYPE inconnu --> force a REAL'
381           status = NF90_DEF_VAR(kcdf_id,ycdfvar,type_float,&
382                    ivdims(:invdims),tpreclist(ji)%id)
383           IF (status /= NF90_NOERR) CALL HANDLE_ERR(status,__LINE__)
384           
385
386        END SELECT
387
388        ! Compress data (costly operation for the CPU)
389        IF (ocompress .AND. invdims>0) THEN
390          status = NF90_DEF_VAR_DEFLATE(kcdf_id,tpreclist(ji)%id,1,1,compress_level)
391          IF (status /= NF90_NOERR) CALL HANDLE_ERR(status,__LINE__)
392        END IF
393
394        ! GRID attribute definition
395        status = NF90_PUT_ATT(kcdf_id,tpreclist(ji)%id,'GRID',tpreclist(ji)%grid)
396        IF (status /= NF90_NOERR) CALL HANDLE_ERR(status,__LINE__)
397
398        ! COMMENT attribute definition
399        status = NF90_PUT_ATT(kcdf_id,tpreclist(ji)%id,'COMMENT',trim(tpreclist(ji)%comment))
400        IF (status /= NF90_NOERR) CALL HANDLE_ERR(status,__LINE__)
401
402     END DO
403     
404     status = NF90_ENDDEF(kcdf_id)
405     IF (status /= NF90_NOERR) CALL HANDLE_ERR(status,__LINE__)
406
407     
408   END SUBROUTINE def_ncdf
409
410   SUBROUTINE fill_ncdf(klu,kcdf_id,tpreclist,knaf,kbuflen,current_level)
411     INTEGER,                      INTENT(IN):: klu
412     INTEGER,                      INTENT(IN):: kcdf_id
413     TYPE(workfield), DIMENSION(:),INTENT(IN):: tpreclist
414     INTEGER,                      INTENT(IN):: knaf
415     INTEGER,                      INTENT(IN):: kbuflen
416     INTEGER, INTENT(IN), OPTIONAL           :: current_level
417
418 #ifdef LOWMEM
419     INTEGER(KIND=8),DIMENSION(:),ALLOCATABLE :: iwork
420 #endif
421     INTEGER                                  :: ji,jj
422     INTEGER,DIMENSION(:),ALLOCATABLE         :: itab
423     REAL   (KIND=8),DIMENSION(:),ALLOCATABLE :: xtab
424     CHARACTER, DIMENSION(:), ALLOCATABLE     :: ytab
425     INTEGER                                  :: status
426     INTEGER                                  :: extent, ndims
427     INTEGER                                  :: ich
428     INTEGER                                  :: level
429     INTEGER(KIND=LFI_INT)                    :: iresp,ilu,ileng,ipos
430     CHARACTER(LEN=4)                         :: suffix
431
432     !
433     ilu = klu
434     !
435
436     IF (present(current_level)) THEN
437       write(suffix,'(I4.4)') current_level
438       level = current_level
439     ElSE
440       suffix=''
441       level = 1
442     END IF
443
444 #if LOWMEM
445     ALLOCATE(iwork(kbuflen))
446 #endif
447     ALLOCATE(itab(kbuflen))
448     ALLOCATE(xtab(kbuflen))
449
450     DO ji=1,knaf
451 #if LOWMEM
452        CALL LFINFO(iresp,ilu,trim(tpreclist(ji)%name)//trim(suffix),ileng,ipos)
453        CALL LFILEC(iresp,ilu,trim(tpreclist(ji)%name)//trim(suffix),iwork,ileng)
454 #endif
455        IF (ASSOCIATED(tpreclist(ji)%dim)) THEN
456           extent = tpreclist(ji)%dim%len
457           ndims = tpreclist(ji)%dim%ndims
458        ELSE
459           extent = 1
460           ndims = 0
461        END IF
462
463        SELECT CASE(tpreclist(ji)%TYPE)
464        CASE (INT,BOOL)
465 #if LOWMEM
466 ***
467 print *,'lowmem: not tested!!!!!' (to be compared with no low mem version)
468          itab(1:extent) = iwork(3+iwork(2):)
469 #else
470          itab(1:extent) = lfiart(ji)%iwtab(3+lfiart(ji)%iwtab(2):)
471 #endif
472 !TODO: works in all cases??? (X, Y, Z dimensions assumed to be ptdimx,y or z)
473          SELECT CASE(ndims)
474          CASE (0)
475            status = NF90_PUT_VAR(kcdf_id,tpreclist(ji)%id,itab(1))
476          CASE (1)
477            status = NF90_PUT_VAR(kcdf_id,tpreclist(ji)%id,itab(1:extent),count=(/extent/))
478          CASE (2)
479            status = NF90_PUT_VAR(kcdf_id,tpreclist(ji)%id,reshape(itab,(/ptdimx%len,ptdimy%len/)), &
480                                  start = (/1,1,level/) )
481          CASE (3)
482            status = NF90_PUT_VAR(kcdf_id,tpreclist(ji)%id,reshape(itab,(/ptdimx%len,ptdimy%len,ptdimz%len/)))
483          CASE DEFAULT
484            print *,'Error: arrays with ',tpreclist(ji)%dim%ndims,' dimensions are not supported'
485          END SELECT
486          
487        CASE (FLOAT)
488 #if LOWMEM
489 ***
490 print *,'lowmem: not tested!!!!!' (to be compared with no low mem version)
491          xtab(1:extent) = TRANSFER(iwork(3+iwork(2):),(/ 0.0_8 /))
492 #else
493          xtab(1:extent) = TRANSFER(lfiart(ji)%iwtab(3+lfiart(ji)%iwtab(2):),(/ 0.0_8 /))
494 #endif
495 !TODO: works in all cases??? (X, Y, Z dimensions assumed to be ptdimx,y or z)
496          SELECT CASE(ndims)
497          CASE (0)
498            status = NF90_PUT_VAR(kcdf_id,tpreclist(ji)%id,xtab(1))
499          CASE (1)
500            status = NF90_PUT_VAR(kcdf_id,tpreclist(ji)%id,xtab(1:extent),count=(/extent/))
501          CASE (2)
502            status = NF90_PUT_VAR(kcdf_id,tpreclist(ji)%id,reshape(xtab,(/ptdimx%len,ptdimy%len/)), &
503                                  start = (/1,1,level/) )
504          CASE (3)
505            status = NF90_PUT_VAR(kcdf_id,tpreclist(ji)%id,reshape(xtab,(/ptdimx%len,ptdimy%len,ptdimz%len/)))
506          CASE DEFAULT
507            print *,'Error: arrays with ',tpreclist(ji)%dim%ndims,' dimensions are not supported'
508          END SELECT
509
510        CASE (TEXT)
511          ALLOCATE(ytab(extent))
512          DO jj=1,extent
513 #if LOWMEM
514 ***
515 print *,'lowmem: not tested!!!!!' (to be compared with no low mem version)
516            ich = iwork(2+iwork(2)+jj)
517 #else
518            ich = lfiart(ji)%iwtab(2+lfiart(ji)%iwtab(2)+jj)
519 #endif
520            ytab(jj) = CHAR(ich)
521          END DO
522          status = NF90_PUT_VAR(kcdf_id,tpreclist(ji)%id,ytab,count=(/extent/))
523          IF (status /= NF90_NOERR) CALL HANDLE_ERR(status,__LINE__)
524          DEALLOCATE(ytab)
525        CASE default
526 #if LOWMEM
527 ***
528 print *,'lowmem: not tested!!!!!' (to be compared with no low mem version)
529          xtab(1:extent) = TRANSFER(iwork(3+iwork(2):),(/ 0.0_8 /))
530 #else         
531          xtab(1:extent) = TRANSFER(lfiart(ji)%iwtab(3+lfiart(ji)%iwtab(2):),(/ 0.0_8 /))
532 #endif
533 !TODO: works in all cases??? (X, Y, Z dimensions assumed to be ptdimx,y or z)
534          SELECT CASE(ndims)
535          CASE (0)
536            status = NF90_PUT_VAR(kcdf_id,tpreclist(ji)%id,xtab(1))
537          CASE (1)
538            status = NF90_PUT_VAR(kcdf_id,tpreclist(ji)%id,xtab(1:extent),count=(/extent/))
539          CASE (2)
540            status = NF90_PUT_VAR(kcdf_id,tpreclist(ji)%id,reshape(xtab,(/ptdimx%len,ptdimy%len/)), &
541                                  start = (/1,1,level/) )
542          CASE (3)
543            status = NF90_PUT_VAR(kcdf_id,tpreclist(ji)%id,reshape(xtab,(/ptdimx%len,ptdimy%len,ptdimz%len/)))
544          CASE DEFAULT
545            print *,'Error: arrays with ',tpreclist(ji)%dim%ndims,' dimensions are not supported'
546          END SELECT
547          IF (status /= NF90_NOERR) CALL HANDLE_ERR(status,__LINE__)
548
549        END SELECT
550
551     END DO
552     DEALLOCATE(itab,xtab)
553 #if LOWMEM
554     DEALLOCATE(iwork)
555 #endif 
556   END SUBROUTINE fill_ncdf
557
558   SUBROUTINE parse_cdf(kcdf_id,tpreclist,kbuflen)
559     INTEGER, INTENT(IN)                    :: kcdf_id
560     TYPE(workfield), DIMENSION(:), POINTER :: tpreclist
561     INTEGER, INTENT(OUT)                   :: kbuflen
562
563
564     INTEGER :: status
565     INTEGER :: nvars, var_id
566     INTEGER :: jdim
567     INTEGER :: sizemax
568     INTEGER :: itype
569     INTEGER, DIMENSION(10) :: idim_id
570     INTEGER :: icomlen,idimlen,idims,idimtmp
571     
572     status = NF90_INQUIRE(kcdf_id, nvariables = nvars)
573     IF (status /= NF90_NOERR) CALL HANDLE_ERR(status,__LINE__)
574     ALLOCATE(tpreclist(nvars))
575
576     sizemax = 0
577
578     CALL init_dimCDF()
579     
580     ! Parcours de toutes les variables et extraction des infos
581     !      - nom de dimension
582     !      - dimension, etendue
583     !      - attributs
584     DO var_id = 1, nvars
585        ! Pour la forme
586        tpreclist(var_id)%id = var_id  
587        
588        ! Nom, type et dimensions de la variable
589        status = NF90_INQUIRE_VARIABLE(kcdf_id, var_id, name = tpreclist(var_id)%name, xtype = itype, ndims = idims, &
590                                       dimids = idim_id)
591        IF (status /= NF90_NOERR) CALL HANDLE_ERR(status,__LINE__)
592        
593        SELECT CASE(itype)
594        CASE(NF90_CHAR)
595           tpreclist(var_id)%TYPE = TEXT
596        CASE(NF90_INT)
597           tpreclist(var_id)%TYPE = INT
598        CASE(NF90_FLOAT,NF90_DOUBLE)
599           tpreclist(var_id)%TYPE = FLOAT
600        CASE default 
601           PRINT *, 'Attention : variable ',TRIM(tpreclist(var_id)&
602                & %name), ' a un TYPE non reconnu par le convertisseur.'
603           PRINT *, '--> TYPE force a REAL(KIND 8) dans LFI !'
604        END SELECT
605       
606        IF (idims == 0) THEN
607           ! variable scalaire
608           NULLIFY(tpreclist(var_id)%dim)
609           idimlen = 1
610        ELSE
611           ! infos sur dimensions
612           idimlen = 1
613           DO jdim=1,idims
614             status = NF90_INQUIRE_DIMENSION(kcdf_id,idim_id(jdim),len = idimtmp)
615             IF (status /= NF90_NOERR) CALL HANDLE_ERR(status,__LINE__)
616             idimlen = idimlen*idimtmp
617           END DO
618           
619           tpreclist(var_id)%dim=>get_dimCDF(idimlen)
620           ! seul le champ 'len' de dimCDF sera utilise par la suite
621        END IF
622        
623        ! GRID et COMMENT attributes
624        status = NF90_GET_ATT(kcdf_id,var_id,'GRID',tpreclist(var_id)%grid)
625        IF (status /= NF90_NOERR) CALL HANDLE_ERR(status,__LINE__)
626
627        status = NF90_INQUIRE_ATTRIBUTE(kcdf_id,var_id,'COMMENT',len = icomlen)
628        IF (status /= NF90_NOERR) CALL HANDLE_ERR(status,__LINE__)
629        
630        ALLOCATE(character(len=icomlen) :: tpreclist(var_id)%comment)
631        status = NF90_GET_ATT(kcdf_id,var_id,'COMMENT',tpreclist(var_id)%comment)
632        IF (status /= NF90_NOERR) CALL HANDLE_ERR(status,__LINE__)
633
634        
635        IF (sizemax < icomlen+idimlen) sizemax = icomlen+idimlen 
636
637     END DO
638     
639     kbuflen = sizemax
640
641   END SUBROUTINE parse_cdf
642
643   SUBROUTINE build_lfi(kcdf_id,klu,tpreclist,kbuflen)
644     INTEGER,                       INTENT(IN) :: kcdf_id 
645     INTEGER,                       INTENT(IN) :: klu
646     TYPE(workfield), DIMENSION(:), INTENT(IN) :: tpreclist
647     INTEGER,                       INTENT(IN) :: kbuflen
648     
649     INTEGER :: status
650     INTEGER :: ivar,jj
651     INTEGER(KIND=8), DIMENSION(:), POINTER  :: iwork
652     INTEGER(KIND=8), DIMENSION(:), POINTER  :: idata
653     REAL(KIND=8), DIMENSION(:), ALLOCATABLE :: xtab
654     INTEGER,      DIMENSION(:), ALLOCATABLE :: itab
655     CHARACTER,    DIMENSION(:), ALLOCATABLE :: ytab
656     CHARACTER(LEN=FM_FIELD_SIZE)            :: yrecfm
657
658     INTEGER :: iartlen, idlen, icomlen
659     INTEGER(KIND=LFI_INT) :: iresp,ilu,iartlen8
660
661     ! Un article LFI est compose de :
662     !   - 1 entier identifiant le numero de grille
663     !   - 1 entier contenant la taille du commentaire
664     !   - le commentaire code en entier 64 bits
665     !   - les donnees proprement dites
666
667     PRINT *,'Taille buffer = ',2+kbuflen
668
669     ALLOCATE(iwork(2+kbuflen))
670     ALLOCATE(itab(2+kbuflen))
671     ALLOCATE(xtab(2+kbuflen))
672
673     DO ivar=1,SIZE(tpreclist)
674        icomlen = LEN(tpreclist(ivar)%comment)
675
676        ! traitement Grille et Commentaire
677        iwork(1) = tpreclist(ivar)%grid
678        iwork(2) = icomlen
679        DO jj=1,iwork(2)
680           iwork(2+jj)=ICHAR(tpreclist(ivar)%comment(jj:jj))
681        END DO
682
683        IF (ASSOCIATED(tpreclist(ivar)%dim)) THEN
684           idlen = tpreclist(ivar)%dim%len
685        ELSE 
686           idlen = 1
687        END IF
688        
689        iartlen = 2+icomlen+idlen
690        idata=>iwork(3+icomlen:iartlen)
691
692
693        SELECT CASE(tpreclist(ivar)%TYPE)
694        CASE(INT,BOOL)
695           status = NF90_GET_VAR(kcdf_id,tpreclist(ivar)%id,itab)
696           IF (status /= NF90_NOERR) CALL HANDLE_ERR(status,__LINE__)
697
698 !          PRINT *,'INT,BOOL --> ',tpreclist(ivar)%name,',len = ',idlen
699           idata(1:idlen) = itab(1:idlen)
700
701        CASE(FLOAT)
702           status = NF90_GET_VAR(kcdf_id,tpreclist(ivar)%id,xtab)
703           IF (status /= NF90_NOERR) CALL HANDLE_ERR(status,__LINE__)
704           
705 !          PRINT *,'FLOAT    --> ',tpreclist(ivar)%name,',len = ',idlen
706           ! La ligne suivante ne pose aucun pb sur Cray alors que sur
707           ! fuji, elle genere une erreur d'execution
708 !          idata(1:idlen) = TRANSFER(xtab(1:idlen),(/ 0_8 /))
709           
710           ! la correction pour Fuji (valable sur CRAY) est :
711           idata(1:idlen) = TRANSFER(xtab,(/ 0_8 /),idlen)
712
713 !          IF (idlen < 10) PRINT *,'xtab = ',xtab(1:idlen)
714
715        CASE(TEXT)
716           ALLOCATE(ytab(idlen))
717           status = NF90_GET_VAR(kcdf_id,tpreclist(ivar)%id,ytab)
718           IF (status /= NF90_NOERR) CALL HANDLE_ERR(status,__LINE__)
719
720 !          PRINT *,'TEXT -->     ',tpreclist(ivar)%name,',len = ',idlen
721
722           DO jj=1,idlen
723              idata(jj) = ICHAR(ytab(jj))
724           END DO
725           
726           DEALLOCATE(ytab)
727
728        CASE default
729           status = NF90_GET_VAR(kcdf_id,tpreclist(ivar)%id,xtab)
730           IF (status /= NF90_NOERR) CALL HANDLE_ERR(status,__LINE__)
731
732           PRINT *,'Default (ERROR) -->',tpreclist(ivar)%name,',len = ',idlen
733           idata(1:idlen) = TRANSFER(xtab,(/ 0_8 /),idlen)
734
735        END SELECT
736        
737        ! Attention restoration des '%' dans le nom des champs LFI
738        yrecfm = str_replace(tpreclist(ivar)%name,'__','%')
739        ! et des '.'
740        yrecfm = str_replace(yrecfm,'--','.')
741        ilu = klu
742        iartlen8 = iartlen
743        CALL LFIECR(iresp,ilu,yrecfm,iwork,iartlen8)
744
745     END DO
746     DEALLOCATE(iwork,itab,xtab)
747
748   END SUBROUTINE build_lfi
749
750   SUBROUTINE OPEN_FILES(hinfile,houtfile,olfi2cdf,olfilist,ohdf5,kcdf_id,klu,knaf)
751     LOGICAL,          INTENT(IN)  :: olfi2cdf, olfilist, ohdf5
752     CHARACTER(LEN=*), INTENT(IN)  :: hinfile
753     CHARACTER(LEN=*), INTENT(IN)  :: houtfile
754     INTEGER         , INTENT(OUT) :: kcdf_id,klu,knaf
755
756     INTEGER                     :: extindex
757     INTEGER(KIND=LFI_INT)       :: ilu,iresp,iverb,inap,inaf
758     INTEGER                     :: status
759     CHARACTER(LEN=4)            :: ypextsrc, ypextdest
760     LOGICAL                     :: fexist
761     INTEGER                     :: omode
762
763     iverb = 0
764     ilu   = 11
765
766     CALL init_sysfield()
767
768     IF (olfi2cdf) THEN 
769        ! Cas LFI -> NetCDF
770        CALL LFIOUV(iresp,ilu,ltrue,hinfile,'OLD',lfalse&
771             & ,lfalse,iverb,inap,inaf)
772
773        IF (olfilist) THEN
774           CALL LFILAF(iresp,ilu,lfalse)
775           CALL LFIFER(iresp,ilu,'KEEP')
776           return
777        end IF
778
779        IF (ohdf5) THEN
780           status = NF90_CREATE(houtfile, IOR(NF90_CLOBBER,NF90_NETCDF4), kcdf_id)
781        ELSE
782           status = NF90_CREATE(houtfile, IOR(NF90_CLOBBER,NF90_64BIT_OFFSET), kcdf_id) 
783        end IF
784        
785        IF (status /= NF90_NOERR) CALL HANDLE_ERR(status,__LINE__)
786
787        status = NF90_SET_FILL(kcdf_id,NF90_NOFILL,omode)
788        IF (status /= NF90_NOERR) CALL HANDLE_ERR(status,__LINE__)
789 !!$       SELECT CASE(omode)
790 !!$       CASE (NF90_FILL)
791 !!$          PRINT *,'Ancien mode : NF90_FILL'
792 !!$       CASE (NF90_NOFILL)
793 !!$          PRINT *,'Ancien mode : NF90_NOFILL'
794 !!$       CASE default
795 !!$          PRINT *, 'Ancien mode : inconnu'
796 !!$       END SELECT
797        
798     ELSE
799        ! Cas NetCDF -> LFI
800        status = NF90_OPEN(hinfile,NF90_NOWRITE,kcdf_id)
801        IF (status /= NF90_NOERR) CALL HANDLE_ERR(status,__LINE__)
802        
803        inap = 100
804        CALL LFIOUV(iresp,ilu,ltrue,houtfile,'NEW'&
805             & ,lfalse,lfalse,iverb,inap,inaf)
806     END IF
807
808     klu  = ilu
809     knaf = inaf
810
811     PRINT *,'--> Fichier converti : ', houtfile
812
813   END SUBROUTINE OPEN_FILES
814   
815   SUBROUTINE CLOSE_FILES(klu,kcdf_id)
816     INTEGER, INTENT(IN) :: klu, kcdf_id
817     
818     INTEGER(KIND=LFI_INT) :: iresp,ilu
819     INTEGER               :: status
820
821     ilu = klu
822     ! close LFI file
823     CALL LFIFER(iresp,ilu,'KEEP')
824
825     ! close NetCDF file
826     status = NF90_CLOSE(kcdf_id)
827     IF (status /= NF90_NOERR) CALL HANDLE_ERR(status,__LINE__)
828     
829   END SUBROUTINE CLOSE_files
830   
831 END MODULE mode_util