source: trunk/LMDZ.TITAN/libf/muphytitan/asciiread.f90 @ 1793

Last change on this file since 1793 was 1793, checked in by jvatant, 7 years ago

Making Titan's hazy again, part I
+ Added the source folder libf/muphytitan which contains

YAMMS ( Titan's microphysical model ) from J. Burgalat

+ Modif. compilation files linked to this change
JVO

File size: 14.0 KB
Line 
1! Copyright Jérémie Burgalat (2010-2015)
2!
3! burgalat.jeremie@gmail.com
4!
5! This software is a computer program whose purpose is to provide configuration
6! file and command line arguments parsing features to Fortran programs.
7!
8! This software is governed by the CeCILL-B license under French law and
9! abiding by the rules of distribution of free software.  You can  use,
10! modify and/ or redistribute the software under the terms of the CeCILL-B
11! license as circulated by CEA, CNRS and INRIA at the following URL
12! "http://www.cecill.info".
13!
14! As a counterpart to the access to the source code and  rights to copy,
15! modify and redistribute granted by the license, users are provided only
16! with a limited warranty  and the software's author,  the holder of the
17! economic rights,  and the successive licensors  have only  limited
18! liability.
19!
20! In this respect, the user's attention is drawn to the risks associated
21! with loading,  using,  modifying and/or developing or reproducing the
22! software by the user in light of its specific status of free software,
23! that may mean  that it is complicated to manipulate,  and  that  also
24! therefore means  that it is reserved for developers  and  experienced
25! professionals having in-depth computer knowledge. Users are therefore
26! encouraged to load and test the software's suitability as regards their
27! requirements in conditions enabling the security of their systems and/or
28! data to be ensured and,  more generally, to use and operate it in the
29! same conditions as regards security.
30!
31! The fact that you are presently reading this means that you have had
32! knowledge of the CeCILL-B license and that you accept its terms.
33
34!! file: asciiread.F90
35!! summary: ASCII data file reader source file
36!! author: burgalat
37!! date: 2013-2015
38MODULE ASCIIREAD
39  !! ASCII data file reader module
40  !!
41  !! This module provides a single generic method that can be used to read 2D/3D
42  !! data array from ASCII file.
43  !!
44  !! ```fortran
45  !! FUNCTION read_data(path,data) RESULT(err)
46  !! ```
47  !!
48  !! Where:
49  !!
50  !! - __path__ is a string with the path of the data file.
51  !! - __data__ is an output __allocatable__ 2D/3D array of real(kind=8) values.
52  !!
53  !! ## Expected Format of the data file
54  !!
55  !! The input file:
56  !!   - must use blank space(s) as value delimiter.
57  !!   - must have a regular number of columns, that is each data line must
58  !!     have the same number of columns.
59  !!   - can contains any number of empty lines and/or comment line (i.e. line
60  !!     where first non-blank character is "#"). All other lines are assumed
61  !!     to be data.
62  !! Moreover, in the case of 3D data, it must use a SINGLE empty line for "depth" block separator. 
63  !!
64  !! Error occured when:
65  !!
66  !! - path does not refer to a existing file (1)
67  !! - Logical unit 666 is not free (1)
68  !! - the file does not have regular data-columns number (5)
69  !! - at least a value cannot be cast in double precision (5)
70  !!
71  !! On error, __data__ array is not allocated.
72  !!
73  !! ## 3D data structure
74  !!
75  !! In the case of 3D data the method assumes that the input files consists in _D_ blocks
76  !! of _R_ lines with _C_ columns. Each block must be separated by a single empty line
77  !! and each columns must separated by one or more blank spaces (no tabulation ALLOWED).
78  !! 
79  !! On success, the shape of the 3D output array will be _data(R,C,D)_.
80  USE, INTRINSIC :: ISO_FORTRAN_ENV, ONLY: IOSTAT_END, IOSTAT_EOR
81  USE STRINGS, ONLY: tokenize,from_string, st_slen
82  USE ERRORS
83  IMPLICIT NONE
84
85  PRIVATE
86  PUBLIC :: noerror,error, error_to_string,aborting
87  PUBLIC :: read_data, OPERATOR(/=), OPERATOR(==)
88
89  !! Global interface to reading methods
90  INTERFACE read_data
91    MODULE PROCEDURE read_data_2d, read_data_3d
92  END INTERFACE
93
94  CONTAINS
95
96
97  FUNCTION read_data_3d(path,data3d,delimiter) RESULT(err)
98    !! Read an input ASCII data file and stores its content in a 3D array
99    !!
100    !! The function reads an ASCII file and saves its values in a real(kind=8) 3D-array.
101    !!
102    !! The input file:
103    !!
104    !! - must have a regular number of columns, that is each data line must have the same number
105    !!   of columns (according to the delimiter used). 
106    !! - must use a SINGLE empty line for "depth" block separator.
107    !! - can contains any number of comment lines (i.e. line where first non-blank character is "#").
108    !!   All other lines (except empty lines) are assumed to be data.
109    !!
110    !! Error occured when:
111    !!    - Path does not refer to a existing file (-11)
112    !!    - No free logical unit available (-12)
113    !!    - The file does not have regular data-columns number (-5)
114    !!    - At least a value cannot be cast in double precision (-10)
115    !!
116    !! The method assumes the input files consists in _D_ block of _R_ lines
117    !! with _C_ columns. Each block must be separated by a single empty line and
118    !! each columns must separated by one or more blank spaces (no tabulation ALLOWED).
119    !!
120    !! On success, the shape of the 3D output array will be _output(R,C,D)_.
121    !! On error, the 3D output array is __not allocated__.
122    !!
123    !! @note
124    !! The function uses the logical unit 666 !
125    CHARACTER(len=*), INTENT(in)                             :: path      !! Path of the input data file
126    REAL(kind=8), INTENT(out), DIMENSION(:,:,:), ALLOCATABLE :: data3d    !! 3D-array with the output values (double precision)
127    CHARACTER(len=*), INTENT(in), OPTIONAL                   :: delimiter !! Optional column delimiter(s)
128    TYPE(error) :: err                                                    !! Error status of the function
129    LOGICAL                                           :: ok
130    INTEGER                                           :: i,lc,tlc
131    INTEGER                                           :: ndr,ndc,ndd
132    INTEGER                                           :: ir,jc,kd
133    REAL(kind=8), DIMENSION(:), ALLOCATABLE           :: tmp
134    CHARACTER(len=5)                                  :: slc
135    CHARACTER(len=15)                                 :: i2s
136    CHARACTER(len=:), ALLOCATABLE                     :: line,lm1,zdelim
137    CHARACTER(len=st_slen), DIMENSION(:), ALLOCATABLE :: wrds
138
139    zdelim = CHAR(9)//CHAR(32)
140    IF (PRESENT(delimiter)) zdelim = delimiter
141    ! Check file status
142    INQUIRE(FILE=TRIM(path),EXIST=ok)
143    IF (.NOT.ok) THEN
144      err = error(trim(path)//": no such file",-1) ; RETURN
145    ENDIF
146    INQUIRE(unit=666,OPENED=ok)
147    IF (ok) THEN
148      err = error("lun 666 is already used",-1) ; RETURN
149    ENDIF
150    ! Open file
151    OPEN(666,FILE=TRIM(path),STATUS='OLD',ACTION='READ')
152
153    ! First pass :
154    ! ------------
155    !   - get size (rows, columns, depth)
156    !   - check size consistendcy
157    !   - check value type
158    lc = 0 ; tlc = 0
159    ndr = -1 ; ndc = -1 ; ndd = 1
160    DO WHILE(readline(666,line))
161      lm1 = line
162      ! Read the line
163      lc = lc + 1 ; WRITE(slc,'(I5)') lc ; slc = ADJUSTL(slc)
164      ! skip comment line
165      IF (INDEX(TRIM(ADJUSTL(line)),"#") == 1) CYCLE
166      ! An empty line: new 2D block
167      IF (LEN_TRIM(line) == 0) THEN
168        ndd = ndd + 1
169        IF (ndr < 0) THEN
170          ndr = tlc
171        ELSEIF (tlc /= ndr) THEN
172          WRITE(i2s,'(I15)') ndd ; i2s = ADJUSTL(i2s)
173          err = error(trim(path)//":Invalid number of lines in block #"//i2s//"(line "//TRIM(slc)//")",-5)
174          RETURN
175        ENDIF
176        tlc = 0
177        CYCLE
178      ENDIF
179      tlc = tlc + 1
180      ! Splits line in words
181      IF (.NOT.tokenize(line,wrds,zdelim,.true.)) THEN
182        ! cannot tokenize
183        err = error(trim(path)//": Cannot parse line "//TRIM(slc),-5)
184        RETURN
185      ELSEIF (.NOT.from_string(wrds,tmp)) THEN
186        ! cannot cast values
187        err = error(trim(path)//": Cannot cast values at line "//TRIM(slc),-5)
188        RETURN
189      ELSEIF (ndc > 0 .AND. ndc /= SIZE(wrds)) THEN
190        ! current number of columns not equal to last one
191        err = error(trim(path)//": Invalid number of columns (line "//TRIM(slc)//")",-5)
192        RETURN
193      ENDIF
194      IF (ndc == -1) ndc = SIZE(wrds)
195    ENDDO
196
197    ! NOTE:
198    ! there is a possible bug if data file ends by an empty line:
199    ! we will have an extra empty bloc !
200    !   current patch: we save the last line of the file and check it:
201    !     - if we have empty line, we reduce ndd by one.
202    IF (LEN_TRIM(lm1) == 0) ndd = ndd-1
203
204    ! Rewind input data file
205    REWIND(666)
206    ! Allocate memory
207    ALLOCATE(data3d(ndr,ndc,ndd))
208    ir = 0 ; kd = 1 ;
209    DO WHILE(readline(666,line))
210      IF (INDEX(TRIM(ADJUSTL(line)),"#") == 1) CYCLE
211      ir = ir + 1
212      ! empty line update block subscripts
213      IF (LEN_TRIM(line) == 0) THEN
214        kd = kd + 1 ; ir = 0 ; CYCLE
215      ENDIF
216      ok = tokenize(line,wrds,zdelim,.true.)
217      DO jc = 1,ndc ; ok = from_string(wrds(jc),data3d(ir,jc,kd)) ; ENDDO
218    ENDDO
219    CLOSE(666)
220  END FUNCTION read_data_3d
221
222  FUNCTION read_data_2d(path,data2d,delimiter) RESULT(err)
223    !! Read an input ASCII data file and stores its content in a 2D array
224    !!
225    !! The function reads an ASCII file and saves its values in a real(kind=8) 2D-array.
226    !!
227    !! The input file:
228    !!
229    !! - can contains any number of empty lines and/or comment line (i.e. line where first
230    !!   non-blank character is "#"). All other lines are assumed to be data.
231    !! - must have a regular number of columns, that is each data line must have the same
232    !!   number of columns.
233    !! - must use blank space(s) as value delimiter.
234    !!
235    !! Error occured when:
236    !!
237    !! - Path does not refer to a existing file (-1)
238    !! - Logical unit 666 is not free (-1)
239    !! - The file does not have regular data-columns number (-5)
240    !! - At least a value cannot be cast in double precision (-5)
241    !!
242    !! On error, the 2D output array is __not allocated__.
243    !! @note
244    !! The function uses the logical unit 666 !
245    CHARACTER(len=*), INTENT(in)                           :: path      !! Path of the input data file
246    REAL(kind=8), INTENT(out), DIMENSION(:,:), ALLOCATABLE :: data2d    !! 2D-array with the output values (double precision)
247    CHARACTER(len=*), INTENT(in), OPTIONAL                 :: delimiter !! Optional column delimiter(s)
248    TYPE(error) :: err                                                  !! Error status of function.
249    LOGICAL                                           :: ok
250    INTEGER                                           :: i,e,vc,lc
251    INTEGER                                           :: nl,nc
252    REAL(kind=8), DIMENSION(:), ALLOCATABLE           :: tmp
253    CHARACTER(len=5)                                  :: slc
254    CHARACTER(len=:), ALLOCATABLE                     :: line,zdelim
255    CHARACTER(len=st_slen), DIMENSION(:), ALLOCATABLE :: wrds
256
257    zdelim = CHAR(9)//CHAR(32)
258    IF (PRESENT(delimiter)) zdelim = delimiter
259    ! Gets array size
260    INQUIRE(FILE=TRIM(path),EXIST=ok)
261    IF (.NOT.ok) THEN
262      err = error(trim(path)//": no such file",-1) ; RETURN
263    ENDIF
264    INQUIRE(unit=666,OPENED=ok)
265    IF (ok) THEN
266      err = error("lun 666 is already used",-1) ; RETURN
267    ENDIF
268    OPEN(666,FILE=TRIM(path),STATUS='OLD',ACTION='READ')
269    vc = 0 ; lc=0 ; ok = .true.
270    ! Read the file twice :)
271    lc=0 ; vc = 0 ; nc=-1
272    ! First pass : get number of row values and checks everything !
273    DO
274      ! Read the line
275      IF (.NOT.readline(666,line)) EXIT
276      lc = lc + 1 ; WRITE(slc,'(I5)') lc ; slc = ADJUSTL(slc)
277      ! skip empty/comment line
278      IF (INDEX(TRIM(ADJUSTL(line)),"#") == 1.OR.LEN_TRIM(line) == 0) CYCLE
279      ! update row counter
280      vc = vc + 1
281      ! Splits line in words
282      IF (.NOT.tokenize(line,wrds,zdelim,.true.)) THEN
283        ! cannot tokenize
284        err = error(trim(path)//": Cannot parse line "//TRIM(slc),-5)
285        RETURN
286      ELSEIF (.NOT.from_string(wrds,tmp)) THEN
287        ! cannot cast values
288        do i=1,size(wrds) ; write(*,*) trim(wrds(i)) ; enddo
289        write(1111,'(a)') line
290        err = error(trim(path)//": Cannot cast values at line "//TRIM(slc),-5)
291        RETURN
292      ELSEIF (nc > 0 .AND. nc /= SIZE(wrds)) THEN
293        ! current number of columns not equal to last one
294        err = error(trim(path)//": Invalid number of columns (line "//TRIM(slc)//")",-5)
295        RETURN
296      ENDIF
297      IF (nc == -1) nc = SIZE(wrds)
298    ENDDO
299    ! Rewind input data file
300    REWIND(666)
301    nl = vc
302    ! allocate memory
303    ALLOCATE(data2d(nl,nc))
304    ! Second pass : saves values :)
305    vc = 0
306    DO WHILE(vc <= nl)
307      ! Reads the line
308      IF (.NOT.readline(666,line)) EXIT
309      ! Check if we have comment or null string
310      IF (INDEX(TRIM(ADJUSTL(line)),"#") == 1.OR.LEN_TRIM(line) == 0) CYCLE
311      vc = vc + 1
312      ok = tokenize(line,wrds,zdelim,.true.)
313      DO i = 1,nc
314        ok = from_string(wrds(i),data2d(vc,i))
315      ENDDO
316    ENDDO
317    CLOSE(666)
318    RETURN
319  END FUNCTION read_data_2d
320
321  FUNCTION readline(lun,line) RESULT(not_eof)
322    !! Read a complete line
323    !!
324    !! Each time, it is called, the function reads a complete of the file opened in __lun__
325    !! logical unit and returns .false. if EOF has been reached, .true. otherwise.
326    !!
327    !! The function is intended to read a file line by line:
328    !!
329    !! ```fortran
330    !! lu = 1
331    !! open(lu,file="path/to/the/file/to/read")
332    !! l=0   ! A line counter
333    !! DO WHILE(readline(lu,line))
334    !!   l = l + 1
335    !!   WRITE(*,'("L",I2.2,": ",(a))') il,line
336    !! ENDDO
337    !! CLOSE(1)
338    !! ```
339    INTEGER, INTENT(in)                        :: lun  !! Logical unit with the opened file to read.
340    CHARACTER(len=:), ALLOCATABLE, INTENT(out) :: line !! Processed line
341    LOGICAL  :: not_eof                                !! .true. if EOF has NOT been reached yet, .false. otherwise
342    CHARACTER(len=50) :: buf
343    INTEGER           :: e,sz
344    not_eof = .true. ; line = ''
345    DO
346      READ(lun,'(a)',ADVANCE="no",SIZE=sz,IOSTAT=e) buf
347      IF (e == IOSTAT_END) THEN
348        not_eof = .false.
349        IF (sz > 0) line=line//buf(1:sz)
350        EXIT
351      ELSE IF (e == IOSTAT_EOR) THEN
352        line = line//buf(1:sz)
353        EXIT
354      ELSE
355        line = line//buf
356      ENDIF
357    ENDDO
358  END FUNCTION readline
359
360END MODULE
Note: See TracBrowser for help on using the repository browser.