Skip to content
Snippets Groups Projects
read_MNHfile.py 11.3 KiB
Newer Older
  • Learn to ignore specific revisions
  • #!/usr/bin/env python3
    # -*- coding: utf-8 -*-
    """
    Created on Mon Feb 22 10:29:13 2021
    
    @author: rodierq
    """
    
    import netCDF4 as nc
    import numpy as np
    
    def read_netcdf(LnameFiles, Dvar_input, path='.', removeHALO=True):
      Dvar = {}
      for i,nameFiles in enumerate(LnameFiles):
        f_nb = 'f' + str(i+1)
        print('Reading file ' + f_nb)
        print(path + nameFiles)
        if '000' in nameFiles[-6:-3]: #time series file (diachronic)
            theFile = nc.Dataset(path + nameFiles,'r')
            if theFile['MASDEV'][0] <= 54:
                read_TIMESfiles_54(theFile, f_nb, Dvar_input, Dvar)
            else : # 55 groups variables
                read_TIMESfiles_55(theFile, f_nb, Dvar_input, Dvar, removeHALO)
            theFile.close()
        else:
            read_BACKUPfile(nameFiles, f_nb, Dvar_input, Dvar, path=path, removeHALO=removeHALO)
      return Dvar #Return the dic of [files][variables]
    
    def read_BACKUPfile(nameF, ifile, Dvar_input, Dvar_output, path='.', removeHALO=True):
      theFile = nc.Dataset(path + nameF,'r')
      Dvar_output[ifile] = {} #initialize dic for each files 
           
      #  Reading date since beginning of the model run
      Dvar_output[ifile]['date'] =  theFile.variables['time'].units
      Dvar_output[ifile]['time'] =  theFile.variables['time'][0]
    
      for var in Dvar_input[ifile]: #For each files
        #  Read variables
        n_dim = theFile.variables[var].ndim
        name_dim = theFile.variables[var].dimensions
     #   print(var + "  n_dim  = " + str(n_dim))
        
        if (n_dim ==0) or (n_dim == 1 and 'time' in name_dim):  #Scalaires or Variable time
          Dvar_output[ifile][var] = theFile.variables[var][0].data
        else:      
            if removeHALO:
              if n_dim == 1:
                Dvar_output[ifile][var] = theFile.variables[var][1:-1] #Variables 1D
              elif n_dim == 2:
                if theFile.variables[var].shape[0] == 1 and 'size' in name_dim[1]: #Variables 2D with the second dimension not a coordinate (--> list of strings : chemicals)
                    Dvar_output[ifile][var] = theFile.variables[var][0,:] #Variables 2D
                elif theFile.variables[var].shape[0] == 1: #Variables 2D with first dim = 0
                    Dvar_output[ifile][var] = theFile.variables[var][0,1:-1] #Variables 2D
                else:
                  Dvar_output[ifile][var] = theFile.variables[var][1:-1,1:-1] #Variables 2D
              elif n_dim == 5: #Variables time, sizeXX, level, nj, ni (ex: chemical budgets in 3D)
                  Dvar_output[ifile][var] = theFile.variables[var][0, :, 1:-1,:1:-1,1:-1]
              elif n_dim == 4 and 'time' in name_dim and ('level' in name_dim or 'level_w' in name_dim): # time,z,y,x
                  if theFile.variables[var].shape[1] == 1: #Variables 4D time,z,y,x with time=0 z=0
                    Dvar_output[ifile][var] = theFile.variables[var][0,0,1:-1,1:-1] #Variables 2D y/x
                  elif theFile.variables[var].shape[2] == 1: #Variable 2D  (0,zz,0,xx)
                    Dvar_output[ifile][var] = theFile.variables[var][0,1:-1,0,1:-1] #Variables 2D z/y
                  elif theFile.variables[var].shape[3] == 1: #Variable 2D  (0,zz,yy,0)
                    Dvar_output[ifile][var] = theFile.variables[var][0,1:-1,1:-1,0] #Variables 2D z/x
                  ## ATTENTION VARIABLE 1D codé en 4D non faite
                  else: #Variable 3D simple
                    Dvar_output[ifile][var] = theFile.variables[var][0,1:-1,1:-1,1:-1] #Variables time + 3D     
              elif n_dim == 4 and 'time' in name_dim and 'level' not in name_dim and 'level_w' not in name_dim: # time,nb_something,y,x
                   Dvar_output[ifile][var] = theFile.variables[var][0,:,1:-1,1:-1] #Variables 2D y/x
              elif n_dim == 3 and 'time' in name_dim: # time, y, x 
                Dvar_output[ifile][var] = theFile.variables[var][0,1:-1,1:-1]
              else :
                Dvar_output[ifile][var] = theFile.variables[var][1:-1,1:-1,1:-1]  #Variables 3D
            else:
              if n_dim == 1:
                Dvar_output[ifile][var] = theFile.variables[var][:] #Variables 1D  
              elif n_dim == 2:
                if theFile.variables[var].shape[0] == 1 and 'size' in name_dim[1]: #Variables 2D with the second dimension not a coordinate (--> list of strings : chemicals)
                  Dvar_output[ifile][var] = theFile.variables[var][0,:] #Variables 2D
                elif theFile.variables[var].shape[0] == 1: #Variables 2D with first dim = 0
                  Dvar_output[ifile][var] = theFile.variables[var][0,:] #Variables 2D
                else:
                  Dvar_output[ifile][var] = theFile.variables[var][:,:] #Variables 2D
              elif n_dim == 5: #Variables time, sizeXX, level, nj, ni (ex: chemical budgets in 3D)
                  Dvar_output[ifile][var] = theFile.variables[var][0,:,:,:,:]
              elif n_dim == 4: # time,z,y,x
                if theFile.variables[var].shape[1] == 1: #Variables 4D time,z,y,x with time=0 z=0
                  Dvar_output[ifile][var] = theFile.variables[var][0,0,:,:] #Variables 2D y/x
                elif theFile.variables[var].shape[2] == 1: #Variable 2D  (0,zz,0,xx)
                  Dvar_output[ifile][var] = theFile.variables[var][0,:,0,:] #Variables 2D z/y
                elif theFile.variables[var].shape[3] == 1: #Variable 2D  (0,zz,yy,0)
                  Dvar_output[ifile][var] = theFile.variables[var][0,:,:,0] #Variables 2D z/x
                ## ATTENTION VARIABLE 1D codé en 4D non faite
                else: #Variable 3D simple
                  Dvar_output[ifile][var] = theFile.variables[var][0,:,:,:] #Variables time + 3D
              elif n_dim ==3 and name_dim in var.dimensions: # time, y, x
                Dvar_output[ifile][var] = theFile.variables[var][0,:,:]
              else:
                Dvar_output[ifile][var] = theFile.variables[var][:,:,:]  #Variables 3D
            #  For all variables except scalars, change Fill_Value to NaN
            Dvar_output[ifile][var]= np.where(Dvar_output[ifile][var] != -99999.0, Dvar_output[ifile][var], np.nan)
            Dvar_output[ifile][var]= np.where(Dvar_output[ifile][var] != 999.0, Dvar_output[ifile][var], np.nan)
    
      theFile.close()
      return Dvar_output #Return the dic of [files][variables]
    
    def read_TIMESfiles_54(theFile, ifile, Dvar_input, Dvar_output):
        Dvar_output[ifile] = {} #initialize dic for each files 
    
        #  Level variable is automatically read without the Halo
        Dvar_output[ifile]['level'] = theFile.variables['level'][1:-1]
        
        #  Time variable is automatically read (time since begging of the run) from the 1st variable of the asked variable to read
        suffix, name_first_var = remove_PROC(Dvar_input[ifile][0])
        try: #  It is possible that there is only one value (one time) in the .000 file, as such time series are not suitable and the following line can't be executed. The time variable is then not written
            increment = theFile.variables[name_first_var+'___DATIM'][1,-1] - theFile.variables[name_first_var+'___DATIM'][0,-1] #-1 is the last entry of the date (current UTC time in seconds)
            length_time = theFile.variables[name_first_var+'___DATIM'].shape[0]
            Dvar_output[ifile]['time'] = np.arange(increment,increment*(length_time+1),increment)
        except:
            pass
        
        for var in Dvar_input[ifile]: #For each files
            suffix, var_name = remove_PROC(var)
            n_dim = theFile.variables[var].ndim
            name_dim = theFile.variables[var].dimensions
    
            #  First, test if the variable is a dimension/coordinate variable
            if (n_dim ==0):  #  Scalaires variable
                 Dvar_output[ifile][var] = theFile.variables[var][0].data
                 pass
            elif n_dim == 1:
                Dvar_output[ifile][var_name] = theFile.variables[var][1:-1]  #  By default, the Halo is always removed because is not in the other variables in any .000 variable
                pass
            elif n_dim == 2:
                Lsize1 = list_size1(n_dim, name_dim)
                if Lsize1 == [True, False]: Dvar_output[ifile][var_name] = theFile.variables[var][0,1:-1] 
                pass
            
            Lsize1 = list_size1(n_dim, name_dim)
            if Lsize1 == [True, False, False, True, True]: Dvar_output[ifile][var_name] = theFile.variables[var][0,:,:,0,0].T # Need to be transposed here
            if Lsize1 == [True, True, False, True, False]: Dvar_output[ifile][var_name] = theFile.variables[var][0,0,:,0,:]
    
        return Dvar_output #Return the dic of [files][variables]
    
    def read_TIMESfiles_55(theFile, ifile, Dvar_input, Dvar_output, removeHALO=False):
        Dvar_output[ifile] = {} #initialize dic for each files 
          
        for var in Dvar_input[ifile]: #For each var
            suffix, var_name = remove_PROC(var)
    
            try: #  NetCDF4 Variables
                n_dim = theFile.variables[var_name].ndim
                #  First, test if the variable is a dimension/coordinate variable
                if (n_dim ==0):  #  Scalaires variable
                    Dvar_output[ifile][var_name] = theFile.variables[var_name][0].data
                else:
                    if(removeHALO):
                        if n_dim == 1:
                            Dvar_output[ifile][var_name] = theFile.variables[var_name][1:-1]
                        elif n_dim == 2:
                            Dvar_output[ifile][var_name] = theFile.variables[var_name][1:-1,1:-1] 
                        else: 
                            raise NameError("Lecture des variables de dimension sup a 2 pas encore implementees pour fichiers .000")
                    else:
                        if n_dim == 1:
                            Dvar_output[ifile][var_name] = theFile.variables[var_name][:]
                        elif n_dim == 2:
                            Dvar_output[ifile][var_name] = theFile.variables[var_name][:,:] 
                        else: 
                            raise NameError("Lecture des variables de dimension sup a 2 pas encore implementees pour fichiers .000")
            except KeyError: # NetCDF4 Group
                if theFile.groups[var_name].type == 'TLES' : #LES type
                    #  Build the true name of the groups.variables :
                    whites = ' '*(17 - len('(cart)') - len(var_name)) #  TODO : a adapter selon le type de la variable cart ou autres
                    Dvar_output[ifile][var] = theFile.groups[var].variables[var + whites + '(cart)'][:,:].T
                elif theFile.groups[var_name].type == 'CART':  #  Budget CART type
                    shapeVar = theFile.groups[var_name].variables[suffix].shape
                    Ltosqueeze=[] #  Build a tuple with the number of the axis which are 0 dimensions to be removed by np.squeeze
                    if shapeVar[0]==1: Ltosqueeze.append(0)
                    if shapeVar[1]==1: Ltosqueeze.append(1)
                    if shapeVar[2]==1: Ltosqueeze.append(2)
                    if shapeVar[3]==1: Ltosqueeze.append(3)
                    Ltosqueeze=tuple(Ltosqueeze)
                    Dvar_output[ifile][var_name] = np.squeeze(theFile.groups[var_name].variables[suffix][:,:,:,:], axis=Ltosqueeze) 
                else:
                    raise NameError("Type de groups variables not implemented in read_MNHfile.py")
        return Dvar_output #Return the dic of [files][variables]
    
    def list_size1(n_dim, named_dim):
        Lsize1 = []
        for i in range(n_dim):
            if 'size1' == named_dim[i]:
                Lsize1.append(True)
            else:
                Lsize1.append(False)
        return Lsize1
        
    def remove_PROC(var):
        if '___PROC' in var:
            var_name = var[:-8]
            suffix = "" # No need of suffix for MNHVERSION < 550 (suffix is for NetCDF4 group
        elif  '___ENDF' in var or '___INIF' in var or '___AVEF' in var:
            var_name = var[:-7]
            suffix = var[-4:]
        else:
            var_name = var
            suffix = ''
        return suffix, var_name