source: trunk/LMDZ.PLUTO/libf/muphypluto/swift_asciiread.F90

! Copyright (c) Jeremie Burgalat (2013-2022)
! Contributor: Jeremie Burgalat (jeremie.burgalat@univ-reims.fr).
!
! This file is part of SWIFT
!
! Permission is hereby granted, free of charge, to any person obtaining a copy of
! this software and associated documentation files (the "Software"), to deal in
! the Software without restriction, including without limitation the rights to
! use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies of
! the Software, and to permit persons to whom the Software is furnished to do so,
! subject to the following conditions:
!
! The above copyright notice and this permission notice shall be included in all
! copies or substantial portions of the Software.
!
! THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
! IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS
! FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR
! COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER
! IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN
! CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.

!! File:    swift_asciiread.F90
!! Summary: ASCII data file reader source file
!! Author:  J. Burgalat
!! Date:    2013-2022

MODULE SWIFT_ASCIIREAD
    !! ASCII data file reader module
    !!
    !! This module provides a single generic method that can be used to read 2D/3D
    !! data array from ASCII file.
    !!
    !! ```
    !! FUNCTION read_data(path,data) RESULT(err)
    !! ```
    !!
    !! Where:
    !!
    !! - __path__ is a string with the path of the data file.
    !! - __data__ is an output __allocatable__ 2D/3D array of real(kind=8) values.
    !!
    !! ## Expected Format of the data file
    !!
    !! The input file:
    !!   - must use blank space(s) as value delimiter.
    !!   - must have a regular number of columns, that is each data line must
    !!     have the same number of columns.
    !!   - can contains any number of empty lines and/or comment line (i.e. line
    !!     where first non-blank character is "#"). All other lines are assumed
    !!     to be data.
    !! Moreover, in the case of 3D data, it must use a SINGLE empty line for "depth" block separator.
    !!
    !! Error occured when:
    !!
    !! - path does not refer to a existing file (1)
    !! - No free logical unit available (1)
    !! - the file does not have regular data-columns number (5)
    !! - at least a value cannot be cast in double precision (5)
    !!
    !! On error, __data__ array is not allocated.
    !!
    !! ## 3D data structure
    !!
    !! In the case of 3D data the method assumes that the input files consists in _D_ blocks
    !! of _R_ lines with _C_ columns. Each block must be separated by a single empty line
    !! and each columns must separated by one or more blank spaces (no tabulation ALLOWED).
    !!
    !! On success, the shape of the 3D output array will be _data(R,C,D)_.
    USE, INTRINSIC :: ISO_FORTRAN_ENV, ONLY: IOSTAT_END, IOSTAT_EOR
    !USE STRING_OP, ONLY: tokenize,from_string, st_slen
    USE SWIFT_STRING_OP
    USE SWIFT_ERRORS
    IMPLICIT NONE
  
    PRIVATE
    PUBLIC :: noerror,error, error_to_string,aborting
    PUBLIC :: readline,read_data, OPERATOR(/=), OPERATOR(==)
  
    !! Global interface to reading methods
    INTERFACE read_data
      MODULE PROCEDURE read_data_2d, read_data_3d
    END INTERFACE
  
    CONTAINS
  
  
    FUNCTION read_data_3d(path,data3d,delimiter) RESULT(err)
      !! Read an input ASCII data file and stores its content in a 3D array
      !!
      !! The function reads an ASCII file and saves its values in a real(kind=8) 3D-array.
      !!
      !! The input file:
      !!
      !! - must have a regular number of columns, that is each data line must have the same number
      !!   of columns (according to the delimiter used).
      !! - must use a SINGLE empty line for "depth" block separator.
      !! - can contains any number of comment lines (i.e. line where first non-blank character is "#").
      !!   All other lines (except empty lines) are assumed to be data.
      !!
      !! Error occured when:
      !!    - Path does not refer to a existing file (-11)
      !!    - No free logical unit available (-1)
      !!    - The file does not have regular data-columns number (-5)
      !!    - At least a value cannot be cast in double precision (-10)
      !!
      !! The method assumes the input files consists in _D_ block of _R_ lines
      !! with _C_ columns. Each block must be separated by a single empty line and
      !! each columns must separated by one or more blank spaces (no tabulation ALLOWED).
      !!
      !! On success, the shape of the 3D output array will be _output(R,C,D)_.
      !! On error, the 3D output array is __not allocated__.
      CHARACTER(len=*), INTENT(in)                             :: path      !! Path of the input data file
      REAL(kind=8), INTENT(out), DIMENSION(:,:,:), ALLOCATABLE :: data3d    !! 3D-array with the output values (double precision)
      CHARACTER(len=*), INTENT(in), OPTIONAL                   :: delimiter !! Optional column delimiter(s)
      TYPE(error) :: err                                                    !! Error status of the function
      LOGICAL                                           :: ok
      INTEGER                                           :: i,lc,tlc
      INTEGER                                           :: ndr,ndc,ndd
      INTEGER                                           :: ir,jc,kd,lu
      REAL(kind=8), DIMENSION(:), ALLOCATABLE           :: tmp
      CHARACTER(len=5)                                  :: slc
      CHARACTER(len=15)                                 :: i2s
      CHARACTER(len=:), ALLOCATABLE                     :: line,lm1,zdelim
      CHARACTER(len=st_slen), DIMENSION(:), ALLOCATABLE :: wrds
  
      zdelim = CHAR(9)//CHAR(32)
      IF (PRESENT(delimiter)) zdelim = delimiter
      ! Check file status
      INQUIRE(FILE=TRIM(path),EXIST=ok)
      IF (.NOT.ok) THEN
        err = error(trim(path)//": no such file",-1) ; RETURN
      ENDIF
      lu = free_lun()
      IF (lu == -1) THEN ; err = error("Cannot find available logical unit...",-1) ; RETURN ; ENDIF
      ! Open file
      OPEN(lu,FILE=TRIM(path),STATUS='OLD',ACTION='READ')
  
      ! First pass :
      ! ------------
      !   - get size (rows, columns, depth)
      !   - check size consistendcy
      !   - check value type
      lc = 0 ; tlc = 0
      ndr = -1 ; ndc = -1 ; ndd = 1
      DO WHILE(readline(lu,line))
        lm1 = line
        ! Read the line
        lc = lc + 1 ; WRITE(slc,'(I5)') lc ; slc = ADJUSTL(slc)
        ! skip comment line
        IF (INDEX(TRIM(ADJUSTL(line)),"#") == 1) CYCLE
        ! An empty line: new 2D block
        IF (LEN_TRIM(line) == 0) THEN
          ndd = ndd + 1
          IF (ndr < 0) THEN
            ndr = tlc
          ELSEIF (tlc /= ndr) THEN
            WRITE(i2s,'(I15)') ndd ; i2s = ADJUSTL(i2s)
            err = error(trim(path)//":Invalid number of lines in block #"//i2s//"(line "//TRIM(slc)//")",-5)
            RETURN
          ENDIF
          tlc = 0
          CYCLE
        ENDIF
        tlc = tlc + 1
        ! Splits line in words
        IF (.NOT.tokenize(line,wrds,zdelim,.true.)) THEN
          ! cannot tokenize
          err = error(trim(path)//": Cannot parse line "//TRIM(slc),-5)
          RETURN
        ELSEIF (.NOT.from_string(wrds,tmp)) THEN
          ! cannot cast values
          err = error(trim(path)//": Cannot cast values at line "//TRIM(slc),-5)
          RETURN
        ELSEIF (ndc > 0 .AND. ndc /= SIZE(wrds)) THEN
          ! current number of columns not equal to last one
          err = error(trim(path)//": Invalid number of columns (line "//TRIM(slc)//")",-5)
          RETURN
        ENDIF
        IF (ndc == -1) ndc = SIZE(wrds)
        IF (ALLOCATED(wrds)) DEALLOCATE(wrds)
        IF (ALLOCATED(tmp)) DEALLOCATE(tmp)
      ENDDO
  
      ! NOTE:
      ! there is a possible bug if data file ends by an empty line:
      ! we will have an extra empty bloc !
      !   current patch: we save the last line of the file and check it:
      !     - if we have empty line, we reduce ndd by one.
      IF (LEN_TRIM(lm1) == 0) ndd = ndd-1
  
      ! Rewind input data file
      REWIND(lu)
      ! Allocate memory
      ALLOCATE(data3d(ndr,ndc,ndd))
      ir = 0 ; kd = 1 ;
      DO WHILE(readline(lu,line))
        IF (INDEX(TRIM(ADJUSTL(line)),"#") == 1) CYCLE
        ir = ir + 1
        ! empty line update block subscripts
        IF (LEN_TRIM(line) == 0) THEN
          kd = kd + 1 ; ir = 0 ; CYCLE
        ENDIF
        ok = tokenize(line,wrds,zdelim,.true.)
        DO jc = 1,ndc ; ok = from_string(wrds(jc),data3d(ir,jc,kd)) ; ENDDO
        IF (ALLOCATED(wrds)) DEALLOCATE(wrds)
      ENDDO
      CLOSE(lu)
    END FUNCTION read_data_3d
  
    FUNCTION read_data_2d(path,data2d,delimiter) RESULT(err)
      !! Read an input ASCII data file and stores its content in a 2D array
      !!
      !! The function reads an ASCII file and saves its values in a real(kind=8) 2D-array.
      !!
      !! The input file:
      !!
      !! - can contains any number of empty lines and/or comment line (i.e. line where first
      !!   non-blank character is "#"). All other lines are assumed to be data.
      !! - must have a regular number of columns, that is each data line must have the same
      !!   number of columns.
      !! - must use blank space(s) as value delimiter.
      !!
      !! Error occured when:
      !!
      !! - Path does not refer to a existing file (-1)
      !! - No free logical unit available (-1)
      !! - The file does not have regular data-columns number (-5)
      !! - At least a value cannot be cast in double precision (-5)
      !!
      !! On error, the 2D output array is __not allocated__.
      USE SWIFT_FSYSTEM
      CHARACTER(len=*), INTENT(in)                           :: path      !! Path of the input data file
      REAL(kind=8), INTENT(out), DIMENSION(:,:), ALLOCATABLE :: data2d    !! 2D-array with the output values (double precision)
      CHARACTER(len=*), INTENT(in), OPTIONAL                 :: delimiter !! Optional column delimiter(s)
      TYPE(error) :: err                                                  !! Error status of function.
      LOGICAL                                           :: ok
      INTEGER                                           :: i,e,vc,lc
      INTEGER                                           :: nl,nc,lu
      REAL(kind=8), DIMENSION(:), ALLOCATABLE           :: tmp
      CHARACTER(len=5)                                  :: slc
      CHARACTER(len=:), ALLOCATABLE                     :: line,zdelim
      CHARACTER(len=st_slen), DIMENSION(:), ALLOCATABLE :: wrds
      zdelim = CHAR(9)//CHAR(32)
      IF (PRESENT(delimiter)) zdelim = delimiter
      INQUIRE(FILE=TRIM(path),EXIST=ok)
      IF (.NOT.ok) THEN
        err = error(trim(path)//": no such file",-1) ; RETURN
      ENDIF
      lu = free_lun()
      IF (lu == -1) THEN ; err = error("Cannot find available logical unit...",-1) ; RETURN ; ENDIF
      OPEN(lu,FILE=TRIM(path),STATUS='OLD',ACTION='READ')
      vc = 0 ; lc=0 ; ok = .true.
      ! Read the file twice :)
      lc=0 ; vc = 0 ; nc=-1
      ! First pass : get number of row values and checks everything !
      DO
        ! Read the line
        IF (.NOT.readline(lu,line)) EXIT
        lc = lc + 1 ; WRITE(slc,'(I5)') lc ; slc = ADJUSTL(slc)
        ! skip empty/comment line
        IF (INDEX(TRIM(ADJUSTL(line)),"#") == 1.OR.LEN_TRIM(line) == 0) CYCLE
        ! update row counter
        vc = vc + 1
        IF (.NOT.tokenize(line,wrds,zdelim,.true.)) THEN
          ! cannot tokenize
          err = error(trim(path)//": Cannot parse line "//TRIM(slc),-5)
          RETURN
        ELSEIF (.NOT.from_string(wrds,tmp)) THEN
          ! cannot cast values
          err = error(trim(path)//": Cannot cast values at line "//TRIM(slc),-5)
          RETURN
        ELSEIF (nc > 0 .AND. nc /= SIZE(tmp)) THEN
          ! current number of columns not equal to last one
          err = error(trim(path)//": Invalid number of columns (line "//TRIM(slc)//")",-5)
          RETURN
        ENDIF
        IF (nc == -1) nc = SIZE(wrds)
        IF (ALLOCATED(wrds)) DEALLOCATE(wrds)
        IF (ALLOCATED(tmp)) DEALLOCATE(tmp)
      ENDDO
      ! Rewind input data file
      REWIND(lu)
      nl = vc
      ! allocate memory
      ALLOCATE(data2d(nl,nc))
      ! Second pass : saves values :)
      vc = 0
      DO WHILE(vc <= nl)
        ! Reads the line
        IF (.NOT.readline(lu,line)) EXIT
        ! Check if we have comment or null string
        IF (INDEX(TRIM(ADJUSTL(line)),"#") == 1.OR.LEN_TRIM(line) == 0) CYCLE
        vc = vc + 1
        ok = tokenize(line,wrds,zdelim,.true.)
        DO i = 1,nc
          ok = from_string(wrds(i),data2d(vc,i))
        ENDDO
        IF (ALLOCATED(wrds)) DEALLOCATE(wrds)
      ENDDO
      CLOSE(lu)
      RETURN
    END FUNCTION read_data_2d
  
    FUNCTION readline(lun,line) RESULT(not_eof)
      !! Read a complete line
      !!
      !! Each time, it is called, the function reads a complete of the file opened in __lun__
      !! logical unit and returns .false. if EOF has been reached, .true. otherwise.
      !!
      !! The function is intended to read a file line by line:
      !!
      !! ```
      !! lu = 1
      !! open(lu,file="path/to/the/file/to/read")
      !! l=0   ! A line counter
      !! DO WHILE(readline(lu,line))
      !!   l = l + 1
      !!   WRITE(*,'("L",I2.2,": ",(a))') il,line
      !! ENDDO
      !! CLOSE(1)
      !! ```
      INTEGER, INTENT(in)                        :: lun  !! Logical unit with the opened file to read.
      CHARACTER(len=:), ALLOCATABLE, INTENT(out) :: line !! Processed line
      LOGICAL  :: not_eof                                !! .true. if EOF has NOT been reached yet, .false. otherwise
      CHARACTER(len=50) :: buf
      INTEGER           :: e,sz
      not_eof = .true. ; line = ''
      DO
        READ(lun,'(a)',ADVANCE="no",SIZE=sz,IOSTAT=e) buf
        IF (e == IOSTAT_END) THEN
          not_eof = .false.
          IF (sz > 0) line=line//buf(1:sz)
          EXIT
        ELSE IF (e == IOSTAT_EOR) THEN
          line = line//buf(1:sz)
          EXIT
        ELSE
          line = line//buf
        ENDIF
      ENDDO
    END FUNCTION readline
  
  END MODULE SWIFT_ASCIIREAD