Changeset 962 for IOIPSL


Ignore:
Timestamp:
03/30/10 11:36:17 (14 years ago)
Author:
bellier
Message:
  • New interface for histbeg (regular 1d or 2D, irregular)
  • New interface for histhori (regular 1d or 2D, irregular)
  • New version of histsync
  • New version of histclo
  • New handling of opened files (using a stack) : the limitation is no more the total number of opened files in a run but is the number of SIMULTANEOUS opened files
File:
1 edited

Legend:

Unmodified
Added
Removed
  • IOIPSL/trunk/src/histcom.f90

    r957 r962  
    3737!- 
    3838  INTERFACE histbeg 
    39     MODULE PROCEDURE histbeg_totreg,histbeg_regular,histbeg_irregular 
     39    MODULE PROCEDURE histb_reg1d,histb_reg2d,histb_irreg 
    4040  END INTERFACE 
    4141!- 
    4242  INTERFACE histhori 
    43     MODULE PROCEDURE histhori_regular,histhori_irregular 
     43    MODULE PROCEDURE histh_reg1d,histh_reg2d,histh_irreg 
    4444  END INTERFACE 
    4545!- 
     
    129129  INTEGER :: tid,xid,yid 
    130130!-General definitions in the NETCDF file 
    131   INTEGER,DIMENSION(2) :: full_size=0,slab_ori,slab_sz 
     131  INTEGER,DIMENSION(2) :: full_size=0,slab_ori,slab_siz 
    132132!-The horizontal axes 
    133133  CHARACTER(LEN=25),DIMENSION(nb_hax_max,2) :: hax_name 
     
    145145TYPE(T_D_F),DIMENSION(nb_files_max),SAVE :: W_F 
    146146!- 
    147 ! Counter of elements 
    148 !- 
    149   INTEGER,SAVE :: nb_files=0 
    150 !- 
    151147! A list of functions which require special action 
    152148! (Needs to be updated when functions are added 
     
    162158!=== 
    163159!- 
    164 SUBROUTINE histbeg_totreg                 & 
    165  & (pfilename,pim,plon,pjm,plat,          & 
    166  &  par_orix,par_szx,par_oriy,par_szy,    & 
    167  &  pitau0,pdate0,pdeltat,phoriid,idf,domain_id,nb_bits) 
    168 !--------------------------------------------------------------------- 
    169 !- This is just an interface for histbeg_regular in case when 
    170 !- the user provides plon and plat as vectors. 
    171 !- Obviously this can only be used for very regular grids. 
    172 !- 
    173 !- INPUT 
    174 !- 
    175 !- pfilename : Name of the netcdf file to be created 
    176 !- pim       : Size of arrays in longitude direction 
    177 !- plon      : Coordinates of points in longitude 
    178 !- pjm       : Size of arrays in latitude direction 
    179 !- plat      : Coordinates of points in latitude 
    180 !- 
    181 !- The next 4 arguments allow to define a horizontal zoom 
    182 !- for this file. It is assumed that all variables to come 
    183 !- have the same index space. This can not be assumed for 
    184 !- the z axis and thus we define the zoom in histdef. 
    185 !- 
    186 !- par_orix  : Origin of the slab of data within the X axis (pim) 
    187 !- par_szx   : Size of the slab of data in X 
    188 !- par_oriy  : Origin of the slab of data within the Y axis (pjm) 
    189 !- par_szy   : Size of the slab of data in Y 
    190 !- 
    191 !- pitau0    : time step at which the history tape starts 
    192 !- pdate0    : The Julian date at which the itau was equal to 0 
    193 !- pdeltat   : Time step, in seconds, of the counter itau 
    194 !-             used in histwrite for instance 
    195 !- 
    196 !- OUTPUT 
    197 !- 
    198 !- phoriid   : ID of the horizontal grid 
    199 !- idf       : ID of the netcdf file 
    200 !- 
    201 !- Optional INPUT arguments 
    202 !- 
    203 !- domain_id       : Domain identifier 
    204 !- 
    205 !- TO DO 
    206 !- 
    207 !- This package should be written in f90 
    208 !- and use the following features : 
    209 !-    - structures for the meta-data of the files and variables 
    210 !-    - memory allocation as needed 
    211 !-    - Pointers 
    212 !- 
    213 !- VERSION 
    214 !- 
     160SUBROUTINE histb_reg1d                 & 
     161 & (pfilename,pim,plon,pjm,plat,       & 
     162 &  par_orix,par_szx,par_oriy,par_szy, & 
     163 &  pitau0,pdate0,pdeltat,phoriid,idf,domain_id,mode) 
     164!--------------------------------------------------------------------- 
     165!- histbeg for 1D regular horizontal coordinates (see histb_all) 
    215166!--------------------------------------------------------------------- 
    216167  IMPLICIT NONE 
     
    224175  REAL,INTENT(IN) :: pdate0,pdeltat 
    225176  INTEGER,INTENT(OUT) :: idf,phoriid 
    226   INTEGER,INTENT(IN),OPTIONAL :: domain_id,nb_bits 
    227 !- 
    228   REAL,ALLOCATABLE,DIMENSION(:,:) :: lon_tmp,lat_tmp 
    229   LOGICAL :: l_dbg 
    230 !--------------------------------------------------------------------- 
    231   CALL ipsldbg (old_status=l_dbg) 
    232 !- 
    233   IF (l_dbg) WRITE(*,*) "histbeg_totreg" 
    234 !- 
    235   ALLOCATE(lon_tmp(pim,pjm),lat_tmp(pim,pjm)) 
    236 !- 
    237   lon_tmp(:,:) = SPREAD(plon(:),2,pjm) 
    238   lat_tmp(:,:) = SPREAD(plat(:),1,pim) 
    239 !- 
    240   CALL histbeg_regular & 
    241  &  (pfilename,pim,lon_tmp,pjm,lat_tmp, & 
    242  &   par_orix,par_szx,par_oriy,par_szy, & 
    243  &   pitau0,pdate0,pdeltat,phoriid,idf, & 
    244  &   .TRUE.,domain_id,nb_bits) 
    245 !- 
    246   DEALLOCATE(lon_tmp,lat_tmp) 
    247 !---------------------------- 
    248 END SUBROUTINE histbeg_totreg 
     177  INTEGER,INTENT(IN),OPTIONAL :: domain_id 
     178  CHARACTER(LEN=*),OPTIONAL,INTENT(IN) :: mode 
     179!--------------------------------------------------------------------- 
     180  CALL histb_all & 
     181 & (1,pfilename,pim,pjm,pitau0,pdate0,pdeltat,phoriid,idf, & 
     182 &  x_1d=plon,y_1d=plat, & 
     183 &  k_orx=par_orix,k_szx=par_szx,k_ory=par_oriy,k_szy=par_szy, & 
     184 &  domain_id=domain_id,mode=mode) 
     185!------------------------- 
     186END SUBROUTINE histb_reg1d 
    249187!=== 
    250 SUBROUTINE histbeg_regular & 
    251  & (pfilename,pim,plon,pjm,plat,           & 
    252  &  par_orix,par_szx,par_oriy,par_szy,     & 
    253  &  pitau0,pdate0,pdeltat,phoriid,idf, & 
    254  &  opt_rectilinear,domain_id,nb_bits) 
    255 !--------------------------------------------------------------------- 
    256 !- This subroutine initializes a netcdf file and returns the ID. 
    257 !- It will set up the geographical space on which the data will be 
    258 !- stored and offers the possibility of seting a zoom. 
    259 !- It also gets the global parameters into the I/O subsystem. 
    260 !- 
    261 !- INPUT 
    262 !- 
    263 !- pfilename : Name of the netcdf file to be created 
    264 !- pim       : Size of arrays in longitude direction 
    265 !- plon      : Coordinates of points in longitude 
    266 !- pjm       : Size of arrays in latitude direction 
    267 !- plat      : Coordinates of points in latitude 
    268 !- 
    269 !- The next 4 arguments allow to define a horizontal zoom 
    270 !- for this file. It is assumed that all variables to come 
    271 !- have the same index space. This can not be assumed for 
    272 !- the z axis and thus we define the zoom in histdef. 
    273 !- 
    274 !- par_orix  : Origin of the slab of data within the X axis (pim) 
    275 !- par_szx   : Size of the slab of data in X 
    276 !- par_oriy  : Origin of the slab of data within the Y axis (pjm) 
    277 !- par_szy   : Size of the slab of data in Y 
    278 !- 
    279 !- pitau0    : time step at which the history tape starts 
    280 !- pdate0    : The Julian date at which the itau was equal to 0 
    281 !- pdeltat   : Time step, in seconds, of the counter itau 
    282 !-             used in histwrite for instance 
    283 !- 
    284 !- OUTPUT 
    285 !- 
    286 !- phoriid   : ID of the horizontal grid 
    287 !- idf       : ID of the netcdf file 
    288 !- 
    289 !- Optional INPUT arguments 
    290 !- 
    291 !- opt_rectilinear : If true we know the grid is rectilinear 
    292 !- domain_id       : Domain identifier 
    293 !- 
    294 !- TO DO 
    295 !- 
    296 !- This package should be written in F90 and use the following 
    297 !- feature : 
    298 !-   - structures for the meta-data of the files and variables 
    299 !-   - memory allocation as needed 
    300 !-   - Pointers 
    301 !- 
    302 !- VERSION 
    303 !- 
     188SUBROUTINE histb_reg2d & 
     189 & (pfilename,pim,plon,pjm,plat,       & 
     190 &  par_orix,par_szx,par_oriy,par_szy, & 
     191 &  pitau0,pdate0,pdeltat,phoriid,idf,domain_id,mode) 
     192!--------------------------------------------------------------------- 
     193!- histbeg for 2D regular horizontal coordinates (see histb_all) 
    304194!--------------------------------------------------------------------- 
    305195  IMPLICIT NONE 
     
    312202  REAL,INTENT(IN) :: pdate0,pdeltat 
    313203  INTEGER,INTENT(OUT) :: idf,phoriid 
    314   LOGICAL,INTENT(IN),OPTIONAL :: opt_rectilinear 
    315   INTEGER,INTENT(IN),OPTIONAL :: domain_id,nb_bits 
    316 !- 
    317   INTEGER :: nfid,iret,m_c 
    318   CHARACTER(LEN=120) :: file 
    319   CHARACTER(LEN=30) :: timenow 
    320   LOGICAL :: rectilinear 
    321   LOGICAL :: l_dbg 
    322 !--------------------------------------------------------------------- 
    323   CALL ipsldbg (old_status=l_dbg) 
    324 !- 
    325   IF (l_dbg) WRITE(*,*) "histbeg_regular 0.0" 
    326 !- 
    327   nb_files = nb_files+1 
    328   IF (nb_files > nb_files_max) THEN 
    329     CALL ipslerr (3,"histbeg", & 
    330    &  'Table of files too small. You should increase nb_files_max', & 
    331    &  'in histcom.f90 in order to accomodate all these files',' ') 
    332   ENDIF 
    333   idf = nb_files 
    334 !- 
    335 ! 1.0 Transfering into the common for future use 
    336 !- 
    337   IF (l_dbg) WRITE(*,*) "histbeg_regular 1.0" 
    338 !- 
    339   W_F(idf)%itau0  = pitau0 
    340   W_F(idf)%date0  = pdate0 
    341   W_F(idf)%deltat = pdeltat 
    342 !- 
    343   IF (PRESENT(opt_rectilinear)) THEN 
    344     rectilinear = opt_rectilinear 
    345   ELSE 
    346     rectilinear = .FALSE. 
    347   ENDIF 
    348 !- 
    349 ! 2.0 Initializes all variables for this file 
    350 !- 
    351   IF (l_dbg) WRITE(*,*) "histbeg_regular 2.0" 
    352 !- 
    353   W_F(idf)%n_var = 0 
    354   W_F(idf)%n_tax = 0 
    355   W_F(idf)%n_hax = 0 
    356   W_F(idf)%n_zax = 0 
    357 !- 
    358   W_F(idf)%slab_ori(1:2) = (/ par_orix,par_oriy /) 
    359   W_F(idf)%slab_sz(1:2)  = (/ par_szx, par_szy  /) 
    360 !- 
    361 ! 3.0 Opening netcdf file and defining dimensions 
    362 !- 
    363   IF (l_dbg) WRITE(*,*) "histbeg_regular 3.0" 
    364 !- 
    365 ! Add DOMAIN number and ".nc" suffix in file name if needed 
    366 !- 
    367   file = pfilename 
    368   CALL flio_dom_file (file,domain_id) 
    369 !- 
    370 ! Check the mode 
    371   IF (PRESENT(nb_bits)) THEN 
    372     SELECT CASE (nb_bits) 
    373     CASE(32) 
    374       m_c = NF90_CLOBBER 
    375     CASE(64) 
    376       m_c = IOR(NF90_CLOBBER,NF90_64BIT_OFFSET) 
    377     CASE DEFAULT 
    378       CALL ipslerr (3,"histbeg", & 
    379  &      'Invalid argument nb_bits for file :',TRIM(file), & 
    380  &      'Supported values are 32 or 64') 
    381     END SELECT 
    382   ELSE 
    383     m_c = IOR(NF90_CLOBBER,NF90_64BIT_OFFSET) 
    384   ENDIF 
    385 !- 
    386 ! Create file 
    387   iret = NF90_CREATE (file,m_c,nfid) 
    388 !- 
    389   IF (rectilinear) THEN 
    390     iret = NF90_DEF_DIM (nfid,'lon',par_szx,W_F(idf)%xid) 
    391     iret = NF90_DEF_DIM (nfid,'lat',par_szy,W_F(idf)%yid) 
    392   ELSE 
    393     iret = NF90_DEF_DIM (nfid,'x',par_szx,W_F(idf)%xid) 
    394     iret = NF90_DEF_DIM (nfid,'y',par_szy,W_F(idf)%yid) 
    395   ENDIF 
    396 !- 
    397 ! 4.0 Declaring the geographical coordinates and other attributes 
    398 !- 
    399   IF (l_dbg) WRITE(*,*) "histbeg_regular 4.0" 
    400 !- 
    401 ! 4.3 Global attributes 
    402 !- 
    403   iret = NF90_PUT_ATT (nfid,NF90_GLOBAL,'Conventions','CF-1.1') 
    404   iret = NF90_PUT_ATT (nfid,NF90_GLOBAL,'file_name',TRIM(file)) 
    405   iret = NF90_PUT_ATT (nfid,NF90_GLOBAL,'production',TRIM(model_name)) 
    406   lock_modname = .TRUE. 
    407   CALL ioget_timestamp (timenow) 
    408   iret = NF90_PUT_ATT (nfid,NF90_GLOBAL,'TimeStamp',TRIM(timenow)) 
    409 !- 
    410 ! 5.0 Saving some important information on this file in the common 
    411 !- 
    412   IF (l_dbg) WRITE(*,*) "histbeg_regular 5.0" 
    413 !- 
    414   IF (PRESENT(domain_id)) THEN 
    415     W_F(idf)%dom_id_svg = domain_id 
    416   ENDIF 
    417   W_F(idf)%ncfid = nfid 
    418   W_F(idf)%full_size(1:2) = (/ pim,pjm /) 
    419 !- 
    420 ! 6.0 storing the geographical coordinates 
    421 !- 
    422   W_F(idf)%regular=.TRUE. 
    423 !- 
    424   CALL histhori_regular (idf,pim,plon,pjm,plat, & 
    425  &  ' ','Default grid',phoriid,rectilinear) 
    426 !----------------------------- 
    427 END SUBROUTINE histbeg_regular 
     204  INTEGER,INTENT(IN),OPTIONAL :: domain_id 
     205  CHARACTER(LEN=*),OPTIONAL,INTENT(IN) :: mode 
     206!--------------------------------------------------------------------- 
     207  CALL histb_all & 
     208 & (2,pfilename,pim,pjm,pitau0,pdate0,pdeltat,phoriid,idf,  & 
     209 &  x_2d=plon,y_2d=plat, & 
     210 &  k_orx=par_orix,k_szx=par_szx,k_ory=par_oriy,k_szy=par_szy,    & 
     211 &  domain_id=domain_id,mode=mode) 
     212!------------------------- 
     213END SUBROUTINE histb_reg2d 
    428214!=== 
    429 SUBROUTINE histbeg_irregular & 
    430  &  (pfilename,pim,plon,plon_bounds,plat,plat_bounds,   & 
    431  &   pitau0,pdate0,pdeltat,phoriid,idf,domain_id,nb_bits) 
    432 !--------------------------------------------------------------------- 
    433 !- This subroutine initializes a netcdf file and returns the ID. 
    434 !- This version is for totaly irregular grids. In this case all 
    435 !- all the data comes in as vectors and for the grid we have 
    436 !- the coordinates of the 4 corners. 
    437 !- It also gets the global parameters into the I/O subsystem. 
    438 !- 
    439 !- INPUT 
    440 !- 
    441 !- pfilename   : Name of the netcdf file to be created 
    442 !- pim         : Size of arrays in longitude direction 
    443 !- plon        : Coordinates of points in longitude 
    444 !- plon_bounds : The 2 corners of the grid in longitude 
    445 !- plat        : Coordinates of points in latitude 
    446 !- plat_bounds : The 2 corners of the grid in latitude 
    447 !- 
    448 !- pitau0    : time step at which the history tape starts 
    449 !- pdate0    : The Julian date at which the itau was equal to 0 
    450 !- pdeltat   : Time step, in seconds, of the counter itau 
    451 !-             used in histwrite for instance 
    452 !- 
    453 !- OUTPUT 
    454 !- 
    455 !- phoriid   : ID of the horizontal grid 
    456 !- idf       : ID of the netcdf file 
    457 !- 
    458 !- Optional INPUT arguments 
    459 !- 
    460 !- domain_id       : Domain identifier 
    461 !- 
    462 !- TO DO 
    463 !- 
    464 !- This package should be written in F90 and use the following 
    465 !- feature : 
    466 !-   - structures for the meta-data of the files and variables 
    467 !-   - memory allocation as needed 
    468 !-   - Pointers 
    469 !- 
    470 !- VERSION 
    471 !- 
     215SUBROUTINE histb_irreg & 
     216 &  (pfilename,pim,plon,plon_bounds,plat,plat_bounds, & 
     217 &   pitau0,pdate0,pdeltat,phoriid,idf,domain_id,mode) 
     218!--------------------------------------------------------------------- 
     219!- histbeg for irregular horizontal coordinates (see histb_all) 
    472220!--------------------------------------------------------------------- 
    473221  IMPLICIT NONE 
     
    480228  REAL,INTENT(IN) :: pdate0,pdeltat 
    481229  INTEGER,INTENT(OUT) :: idf,phoriid 
    482   INTEGER,INTENT(IN),OPTIONAL :: domain_id,nb_bits 
     230  INTEGER,INTENT(IN),OPTIONAL :: domain_id 
     231  CHARACTER(LEN=*),OPTIONAL,INTENT(IN) :: mode 
     232!--------------------------------------------------------------------- 
     233  CALL histb_all & 
     234 & (3,pfilename,pim,pim,pitau0,pdate0,pdeltat,phoriid,idf,  & 
     235 &  x_1d=plon,y_1d=plat,x_bnds=plon_bounds,y_bnds=plat_bounds, & 
     236 &  domain_id=domain_id,mode=mode) 
     237!------------------------- 
     238END SUBROUTINE histb_irreg 
     239!=== 
     240SUBROUTINE histb_all & 
     241 & (k_typ,nc_name,pim,pjm,pitau0,pdate0,pdeltat,phoriid,idf, & 
     242 &  x_1d,y_1d,x_2d,y_2d,k_orx,k_szx,k_ory,k_szy, & 
     243 &  x_bnds,y_bnds,domain_id,mode) 
     244!--------------------------------------------------------------------- 
     245!- General interface for horizontal grids. 
     246!- This subroutine initializes a netcdf file and returns the ID. 
     247!- It will set up the geographical space on which the data will be 
     248!- stored and offers the possibility of seting a zoom. 
     249!- In the case of irregular grids, all the data comes in as vectors 
     250!- and for the grid we have the coordinates of the 4 corners. 
     251!- It also gets the global parameters into the I/O subsystem. 
     252!- 
     253!- INPUT 
     254!- 
     255!- k_typ    : Type of the grid (1 rectilinear, 2 regular, 3 irregular) 
     256!- nc_name  : Name of the netcdf file to be created 
     257!- pim      : Size of arrays in longitude direction 
     258!- pjm      : Size of arrays in latitude direction (pjm=pim for type 3) 
     259!- 
     260!- pitau0   : time step at which the history tape starts 
     261!- pdate0   : The Julian date at which the itau was equal to 0 
     262!- pdeltat  : Time step, in seconds, of the counter itau 
     263!-            used in histwrite for instance 
     264!- 
     265!- OUTPUT 
     266!- 
     267!- phoriid  : Identifier of the horizontal grid 
     268!- idf      : Identifier of the file 
     269!- 
     270!- Optional INPUT arguments 
     271!- 
     272!- For rectilinear or irregular grid 
     273!- x_1d  : The longitudes 
     274!- y_1d  : The latitudes 
     275!- For regular grid 
     276!- x_2d  : The longitudes 
     277!- y_2d  : The latitudes 
     278!- One pair (x_1d,y_1d) or (x_2d,y_2d) must be supplied. 
     279!- 
     280!- For regular grid (reg1d or reg2d),                    
     281!- the next 4 arguments allow to define a horizontal zoom 
     282!- for this file. It is assumed that all variables to come 
     283!- have the same index space. This can not be assumed for 
     284!- the z axis and thus we define the zoom in histdef. 
     285!- k_orx  : Origin of the slab of data within the X axis (pim) 
     286!- k_szx  : Size of the slab of data in X 
     287!- k_ory  : Origin of the slab of data within the Y axis (pjm) 
     288!- k_szy  : Size of the slab of data in Y 
     289!- 
     290!- For irregular grid. 
     291!- x_bnds : The boundaries of the grid in longitude 
     292!- y_bnds : The boundaries of the grid in latitude 
     293!- 
     294!- For all grids. 
     295!- 
     296!- domain_id  : Domain identifier 
     297!- 
     298!- mode       : String of (case insensitive) blank-separated words 
     299!-              defining the mode used to create the file. 
     300!-              Supported keywords : 32, 64 
     301!-              "32/64" defines the offset mode. 
     302!-              The default offset mode is 64 bits. 
     303!-              Keywords "NETCDF4" and "CLASSIC" are reserved 
     304!-              for future use. 
     305!--------------------------------------------------------------------- 
     306  IMPLICIT NONE 
     307!- 
     308  INTEGER,INTENT(IN) :: k_typ 
     309  CHARACTER(LEN=*),INTENT(IN) :: nc_name 
     310  INTEGER,INTENT(IN) :: pim,pjm 
     311  INTEGER,INTENT(IN) :: pitau0 
     312  REAL,INTENT(IN) :: pdate0,pdeltat 
     313  INTEGER,INTENT(OUT) :: idf,phoriid 
     314  REAL,DIMENSION(:),INTENT(IN),OPTIONAL :: x_1d,y_1d 
     315  REAL,DIMENSION(:,:),INTENT(IN),OPTIONAL :: x_2d,y_2d 
     316  INTEGER,INTENT(IN),OPTIONAL :: k_orx,k_szx,k_ory,k_szy 
     317  REAL,DIMENSION(:,:),INTENT(IN),OPTIONAL :: x_bnds,y_bnds 
     318  INTEGER,INTENT(IN),OPTIONAL :: domain_id 
     319  CHARACTER(LEN=*),OPTIONAL,INTENT(IN) :: mode 
    483320!- 
    484321  INTEGER :: nfid,iret,m_c 
    485322  CHARACTER(LEN=120) :: file 
    486323  CHARACTER(LEN=30) :: timenow 
     324  CHARACTER(LEN=11) :: c_nam 
    487325  LOGICAL :: l_dbg 
    488326!--------------------------------------------------------------------- 
    489327  CALL ipsldbg (old_status=l_dbg) 
    490328!- 
    491   IF (l_dbg) WRITE(*,*) "histbeg_irregular 0.0" 
    492 !- 
    493   nb_files = nb_files+1 
    494   IF (nb_files > nb_files_max) THEN 
     329  IF     (k_typ == 1) THEN 
     330    c_nam = 'histb_reg1d' 
     331  ELSEIF (k_typ == 2) THEN 
     332    c_nam = 'histb_reg2d' 
     333  ELSEIF (k_typ == 3) THEN 
     334    c_nam = 'histb_irreg' 
     335  ELSE 
    495336    CALL ipslerr (3,"histbeg", & 
    496  &    'Table of files too small. You should increase nb_files_max', & 
    497  &    'in histcom.f90 in order to accomodate all these files',' ') 
    498   ENDIF 
    499   idf = nb_files 
     337 &    'Illegal value of k_typ argument','in internal interface','?') 
     338  ENDIF 
     339!- 
     340  IF (l_dbg) WRITE(*,*) c_nam//" 0.0" 
     341!- 
     342! Search for a free index 
     343!- 
     344  idf = -1 
     345  DO nfid=1,nb_files_max 
     346    IF (W_F(nfid)%ncfid < 0) THEN 
     347      idf = nfid; EXIT; 
     348    ENDIF 
     349  ENDDO 
     350  IF (idf < 0) THEN 
     351    CALL ipslerr (3,"histbeg", & 
     352   &  'Table of files too small. You should increase nb_files_max', & 
     353   &  'in histcom.f90 in order to accomodate all these files',' ') 
     354  ENDIF 
    500355!- 
    501356! 1.0 Transfering into the common for future use 
    502357!- 
    503   IF (l_dbg) WRITE(*,*) "histbeg_irregular 1.0" 
     358  IF (l_dbg) WRITE(*,*) c_nam//" 1.0" 
    504359!- 
    505360  W_F(idf)%itau0  = pitau0 
     
    509364! 2.0 Initializes all variables for this file 
    510365!- 
    511   IF (l_dbg) WRITE(*,*) "histbeg_irregular 2.0" 
     366  IF (l_dbg) WRITE(*,*) c_nam//" 2.0" 
    512367!- 
    513368  W_F(idf)%n_var = 0 
     
    516371  W_F(idf)%n_zax = 0 
    517372!- 
    518   W_F(idf)%slab_ori(1:2) = (/ 1,1 /) 
    519   W_F(idf)%slab_sz(1:2)  = (/ pim,1 /) 
     373  IF ( (k_typ == 1).OR.(k_typ == 2) ) THEN 
     374    W_F(idf)%slab_ori(1:2) = (/ k_orx,k_ory /) 
     375    W_F(idf)%slab_siz(1:2)  = (/ k_szx,k_szy /) 
     376  ELSE 
     377    W_F(idf)%slab_ori(1:2) = (/ 1,1 /) 
     378    W_F(idf)%slab_siz(1:2) = (/ pim,1 /) 
     379  ENDIF 
    520380!- 
    521381! 3.0 Opening netcdf file and defining dimensions 
    522382!- 
    523   IF (l_dbg) WRITE(*,*) "histbeg_irregular 3.0" 
     383  IF (l_dbg) WRITE(*,*) c_nam//" 3.0" 
    524384!- 
    525385! Add DOMAIN number and ".nc" suffix in file name if needed 
    526386!- 
    527   file  = pfilename 
     387  file = nc_name 
    528388  CALL flio_dom_file (file,domain_id) 
    529389!- 
    530390! Check the mode 
    531   IF (PRESENT(nb_bits)) THEN 
    532     SELECT CASE (nb_bits) 
    533     CASE(32) 
     391!? See fliocom for HDF4 ???????????????????????????????????????????????? 
     392!- 
     393  IF (PRESENT(mode)) THEN 
     394    SELECT CASE (TRIM(mode)) 
     395    CASE('32') 
    534396      m_c = NF90_CLOBBER 
    535     CASE(64) 
     397    CASE('64') 
    536398      m_c = IOR(NF90_CLOBBER,NF90_64BIT_OFFSET) 
    537399    CASE DEFAULT 
    538400      CALL ipslerr (3,"histbeg", & 
    539  &      'Invalid argument nb_bits for file :',TRIM(file), & 
     401 &      'Invalid argument mode for file :',TRIM(file), & 
    540402 &      'Supported values are 32 or 64') 
    541403    END SELECT 
     
    545407!- 
    546408! Create file 
    547   iret = NF90_CREATE (file,m_c,nfid) 
    548 !- 
    549   iret = NF90_DEF_DIM (nfid,'x',pim,W_F(idf)%xid) 
    550   W_F(idf)%yid = 0 
     409!- 
     410  iret = NF90_CREATE(file,m_c,nfid) 
     411!- 
     412  IF     (k_typ == 1) THEN 
     413    iret = NF90_DEF_DIM(nfid,'lon',k_szx,W_F(idf)%xid) 
     414    iret = NF90_DEF_DIM(nfid,'lat',k_szy,W_F(idf)%yid) 
     415  ELSEIF (k_typ == 2) THEN 
     416    iret = NF90_DEF_DIM(nfid,'x',k_szx,W_F(idf)%xid) 
     417    iret = NF90_DEF_DIM(nfid,'y',k_szy,W_F(idf)%yid) 
     418  ELSEIF (k_typ == 3) THEN 
     419    iret = NF90_DEF_DIM(nfid,'x',pim,W_F(idf)%xid) 
     420    W_F(idf)%yid = W_F(idf)%xid 
     421  ENDIF 
    551422!- 
    552423! 4.0 Declaring the geographical coordinates and other attributes 
    553424!- 
    554   IF (l_dbg) WRITE(*,*) "histbeg_irregular 4.0" 
    555 !- 
    556 ! 4.3 Global attributes 
    557 !- 
    558   iret = NF90_PUT_ATT (nfid,NF90_GLOBAL,'Conventions','CF-1.1') 
    559   iret = NF90_PUT_ATT (nfid,NF90_GLOBAL,'file_name',TRIM(file)) 
    560   iret = NF90_PUT_ATT (nfid,NF90_GLOBAL,'production',TRIM(model_name)) 
     425  IF (l_dbg) WRITE(*,*) c_nam//" 4.0" 
     426!- 
     427  iret = NF90_PUT_ATT(nfid,NF90_GLOBAL,'Conventions','CF-1.1') 
     428  iret = NF90_PUT_ATT(nfid,NF90_GLOBAL,'file_name',TRIM(file)) 
     429  iret = NF90_PUT_ATT(nfid,NF90_GLOBAL,'production',TRIM(model_name)) 
    561430  lock_modname = .TRUE. 
    562431  CALL ioget_timestamp (timenow) 
    563   iret = NF90_PUT_ATT (nfid,NF90_GLOBAL,'TimeStamp',TRIM(timenow)) 
     432  iret = NF90_PUT_ATT(nfid,NF90_GLOBAL,'TimeStamp',TRIM(timenow)) 
    564433!- 
    565434! 5.0 Saving some important information on this file in the common 
    566435!- 
    567   IF (l_dbg) WRITE(*,*) "histbeg_irregular 5.0" 
     436  IF (l_dbg) WRITE(*,*) c_nam//" 5.0" 
    568437!- 
    569438  IF (PRESENT(domain_id)) THEN 
     
    571440  ENDIF 
    572441  W_F(idf)%ncfid = nfid 
    573   W_F(idf)%full_size(1:2) = (/ pim,1 /) 
     442  IF ( (k_typ == 1).OR.(k_typ == 2) ) THEN 
     443    W_F(idf)%full_size(1:2) = (/ pim,pjm /) 
     444    W_F(idf)%regular=.TRUE. 
     445  ELSEIF (k_typ == 3) THEN 
     446    W_F(idf)%full_size(1:2) = (/ pim,1 /) 
     447    W_F(idf)%regular=.FALSE. 
     448  ENDIF 
    574449!- 
    575450! 6.0 storing the geographical coordinates 
    576451!- 
    577   W_F(idf)%regular=.FALSE. 
    578 !- 
    579   CALL histhori_irregular & 
    580  &  (idf,pim,plon,plon_bounds,plat,plat_bounds, & 
    581  &   ' ','Default grid',phoriid) 
    582 !------------------------------- 
    583 END SUBROUTINE histbeg_irregular 
     452  IF     (k_typ == 1) THEN 
     453    CALL histh_all & 
     454 &   (k_typ,idf,pim,pjm,' ','Default grid',phoriid, & 
     455 &    x_1d=x_1d,y_1d=y_1d) 
     456  ELSEIF (k_typ == 2) THEN 
     457    CALL histh_all & 
     458 &   (k_typ,idf,pim,pjm,' ','Default grid',phoriid, & 
     459 &    x_2d=x_2d,y_2d=y_2d) 
     460  ELSEIF (k_typ == 3) THEN 
     461    CALL histh_all & 
     462 &   (k_typ,idf,pim,pim,' ','Default grid',phoriid, & 
     463 &    x_1d=x_1d,y_1d=y_1d,x_bnds=x_bnds,y_bnds=y_bnds) 
     464  ENDIF 
     465!----------------------- 
     466END SUBROUTINE histb_all 
    584467!=== 
    585 SUBROUTINE histhori_regular & 
    586  &  (idf,pim,plon,pjm,plat,phname,phtitle,phid,opt_rectilinear) 
    587 !--------------------------------------------------------------------- 
    588 !- This subroutine is made to declare a new horizontale grid. 
     468SUBROUTINE histh_reg1d & 
     469 &  (idf,pim,plon,pjm,plat,phname,phtitle,phid) 
     470!--------------------------------------------------------------------- 
     471!- histhori for 1d regular grid (see histh_all) 
     472!--------------------------------------------------------------------- 
     473  IMPLICIT NONE 
     474!- 
     475  INTEGER,INTENT(IN) :: idf,pim,pjm 
     476  REAL,INTENT(IN),DIMENSION(:) :: plon,plat 
     477  CHARACTER(LEN=*),INTENT(IN) :: phname,phtitle 
     478  INTEGER,INTENT(OUT) :: phid 
     479!--------------------------------------------------------------------- 
     480  CALL histh_all & 
     481 & (1,idf,pim,pjm,phname,phtitle,phid,x_1d=plon,y_1d=plat) 
     482!------------------------- 
     483END SUBROUTINE histh_reg1d 
     484!=== 
     485SUBROUTINE histh_reg2d & 
     486 & (idf,pim,plon,pjm,plat,phname,phtitle,phid) 
     487!--------------------------------------------------------------------- 
     488!- histhori for 2d regular grid (see histh_all) 
     489!--------------------------------------------------------------------- 
     490  IMPLICIT NONE 
     491!- 
     492  INTEGER,INTENT(IN) :: idf,pim,pjm 
     493  REAL,INTENT(IN),DIMENSION(:,:) :: plon,plat 
     494  CHARACTER(LEN=*),INTENT(IN) :: phname,phtitle 
     495  INTEGER,INTENT(OUT) :: phid 
     496!--------------------------------------------------------------------- 
     497  CALL histh_all & 
     498 & (2,idf,pim,pjm,phname,phtitle,phid,x_2d=plon,y_2d=plat) 
     499!------------------------- 
     500END SUBROUTINE histh_reg2d 
     501!=== 
     502SUBROUTINE histh_irreg & 
     503 & (idf,pim,plon,plon_bounds,plat,plat_bounds,phname,phtitle,phid) 
     504!--------------------------------------------------------------------- 
     505!- histhori for irregular grid (see histh_all) 
     506!--------------------------------------------------------------------- 
     507  IMPLICIT NONE 
     508!- 
     509  INTEGER,INTENT(IN) :: idf,pim 
     510  REAL,DIMENSION(:),INTENT(IN) :: plon,plat 
     511  REAL,DIMENSION(:,:),INTENT(IN) :: plon_bounds,plat_bounds 
     512  CHARACTER(LEN=*),INTENT(IN) :: phname,phtitle 
     513  INTEGER,INTENT(OUT) :: phid 
     514!--------------------------------------------------------------------- 
     515  CALL histh_all & 
     516 & (3,idf,pim,pim,phname,phtitle,phid, & 
     517 &  x_1d=plon,y_1d=plat,x_bnds=plon_bounds,y_bnds=plat_bounds) 
     518!------------------------- 
     519END SUBROUTINE histh_irreg 
     520!=== 
     521SUBROUTINE histh_all & 
     522 & (k_typ,idf,pim,pjm,phname,phtitle,phid, & 
     523 &  x_1d,y_1d,x_2d,y_2d,x_bnds,y_bnds) 
     524!--------------------------------------------------------------------- 
     525!- General interface for horizontal grids. 
     526!- This subroutine is made to declare a new horizontal grid. 
    589527!- It has to have the same number of points as 
    590528!- the original and thus in this routine we will only 
     
    596534!- INPUT 
    597535!- 
     536!- k_typ   : Type of the grid (1 rectilinear, 2 regular, 3 irregular) 
    598537!- idf     : The id of the file to which the grid should be added 
    599538!- pim     : Size in the longitude direction 
    600 !- plon    : The longitudes 
    601 !- pjm     : Size in the latitude direction 
    602 !- plat    : The latitudes 
     539!- pjm     : Size in the latitude direction (pjm=pim for type 3) 
    603540!- phname  : The name of grid 
    604541!- phtitle : The title of the grid 
     
    608545!- phid    : Id of the created grid 
    609546!- 
    610 !- OPTIONAL 
    611 !- 
    612 !- opt_rectilinear : If true we know the grid is rectilinear. 
    613 !- 
     547!- Optional INPUT arguments 
     548!- 
     549!- For rectilinear or irregular grid 
     550!- x_1d  : The longitudes 
     551!- y_1d  : The latitudes 
     552!- For regular grid 
     553!- x_2d  : The longitudes 
     554!- y_2d  : The latitudes 
     555!- One pair (x_1d,y_1d) or (x_2d,y_2d) must be supplied. 
     556!- 
     557!- For irregular grid. 
     558!- x_bnds : The boundaries of the grid in longitude 
     559!- y_bnds : The boundaries of the grid in latitude 
    614560!--------------------------------------------------------------------- 
    615561  IMPLICIT NONE 
    616562!- 
     563  INTEGER,INTENT(IN) :: k_typ 
    617564  INTEGER,INTENT(IN) :: idf,pim,pjm 
    618   REAL,INTENT(IN),DIMENSION(pim,pjm) :: plon,plat 
    619565  CHARACTER(LEN=*),INTENT(IN) :: phname,phtitle 
    620566  INTEGER,INTENT(OUT) :: phid 
    621   LOGICAL,INTENT(IN),OPTIONAL :: opt_rectilinear 
     567  REAL,DIMENSION(:),INTENT(IN),OPTIONAL :: x_1d,y_1d 
     568  REAL,DIMENSION(:,:),INTENT(IN),OPTIONAL :: x_2d,y_2d 
     569  REAL,DIMENSION(:,:),INTENT(IN),OPTIONAL :: x_bnds,y_bnds 
    622570!- 
    623571  CHARACTER(LEN=25) :: lon_name,lat_name 
    624   CHARACTER(LEN=80) :: tmp_title,tmp_name 
    625   INTEGER :: ndim 
    626   INTEGER,DIMENSION(2) :: dims 
     572  CHARACTER(LEN=30) :: lonbound_name,latbound_name 
     573  INTEGER :: i_s,i_e 
     574  INTEGER,DIMENSION(2) :: dims,dims_b 
     575  INTEGER :: nbbounds 
     576  INTEGER :: nlonidb,nlatidb,twoid 
     577  LOGICAL :: transp = .FALSE. 
     578  REAL,ALLOCATABLE,DIMENSION(:,:) :: bounds_trans 
     579  REAL :: wmn,wmx 
    627580  INTEGER :: nlonid,nlatid 
    628   INTEGER :: orix,oriy,par_szx,par_szy 
     581  INTEGER :: o_x,o_y,s_x,s_y 
    629582  INTEGER :: iret,nfid 
    630   LOGICAL :: rectilinear 
     583  CHARACTER(LEN=11) :: c_nam 
    631584  LOGICAL :: l_dbg 
    632585!--------------------------------------------------------------------- 
    633586  CALL ipsldbg (old_status=l_dbg) 
    634587!- 
     588  IF     (k_typ == 1) THEN 
     589    c_nam = 'histh_reg1d' 
     590  ELSEIF (k_typ == 2) THEN 
     591    c_nam = 'histh_reg2d' 
     592  ELSEIF (k_typ == 3) THEN 
     593    c_nam = 'histh_irreg' 
     594  ELSE 
     595    CALL ipslerr (3,"histhori", & 
     596 &    'Illegal value of k_typ argument','in internal interface','?') 
     597  ENDIF 
     598!- 
    635599! 1.0 Check that all fits in the buffers 
    636600!- 
    637601  IF (    (pim /= W_F(idf)%full_size(1)) & 
    638  &    .OR.(pjm /= W_F(idf)%full_size(2)) ) THEN 
     602 &    .OR.(W_F(idf)%regular.AND.(pjm /= W_F(idf)%full_size(2)))  & 
     603 &    .OR.(.NOT.W_F(idf)%regular.AND.(W_F(idf)%full_size(2) /= 1)) ) THEN 
    639604    CALL ipslerr (3,"histhori", & 
    640605 &   'The new horizontal grid does not have the same size', & 
     
    643608  ENDIF 
    644609!- 
    645   IF (PRESENT(opt_rectilinear)) THEN 
    646     rectilinear = opt_rectilinear 
    647   ELSE 
    648     rectilinear = .FALSE. 
    649   ENDIF 
    650 !- 
    651610! 1.1 Create all the variables needed 
    652611!- 
    653   IF (l_dbg) WRITE(*,*) "histhori_regular 1.0" 
     612  IF (l_dbg) WRITE(*,*) c_nam//" 1.0" 
    654613!- 
    655614  nfid = W_F(idf)%ncfid 
    656615!- 
    657   ndim = 2 
     616  IF (k_typ == 3) THEN 
     617    IF     (SIZE(x_bnds,DIM=1) == pim) THEN 
     618      nbbounds = SIZE(x_bnds,DIM=2) 
     619      transp = .TRUE. 
     620    ELSEIF (SIZE(x_bnds,DIM=2) == pim) THEN 
     621      nbbounds = SIZE(x_bnds,DIM=1) 
     622      transp = .FALSE. 
     623    ELSE 
     624      CALL ipslerr (3,"histhori", & 
     625 &     'The boundary variable does not have any axis corresponding', & 
     626 &     'to the size of the longitude or latitude variable','.') 
     627    ENDIF 
     628    ALLOCATE(bounds_trans(nbbounds,pim)) 
     629    iret = NF90_DEF_DIM(nfid,'nbnd',nbbounds,twoid) 
     630    dims_b(1:2) = (/ twoid,W_F(idf)%xid /) 
     631  ENDIF 
     632!- 
    658633  dims(1:2) = (/ W_F(idf)%xid,W_F(idf)%yid /) 
    659634!- 
    660   tmp_name = phname 
    661   IF (rectilinear) THEN 
     635  IF     (k_typ == 1) THEN 
    662636    IF (W_F(idf)%n_hax == 0) THEN 
    663637      lon_name = 'lon' 
    664638      lat_name = 'lat' 
    665639    ELSE 
    666       lon_name = 'lon_'//TRIM(tmp_name) 
    667       lat_name = 'lat_'//TRIM(tmp_name) 
    668     ENDIF 
    669   ELSE 
     640      lon_name = 'lon_'//TRIM(phname) 
     641      lat_name = 'lat_'//TRIM(phname) 
     642    ENDIF 
     643  ELSEIF (k_typ == 2) THEN 
    670644    IF (W_F(idf)%n_hax == 0) THEN 
    671645      lon_name = 'nav_lon' 
    672646      lat_name = 'nav_lat' 
    673647    ELSE 
    674       lon_name = 'nav_lon_'//TRIM(tmp_name) 
    675       lat_name = 'nav_lat_'//TRIM(tmp_name) 
    676     ENDIF 
     648      lon_name = 'nav_lon_'//TRIM(phname) 
     649      lat_name = 'nav_lat_'//TRIM(phname) 
     650    ENDIF 
     651  ELSEIF (k_typ == 3) THEN 
     652    IF (W_F(idf)%n_hax == 0) THEN 
     653      lon_name = 'nav_lon' 
     654      lat_name = 'nav_lat' 
     655    ELSE 
     656      lon_name = 'nav_lon_'//TRIM(phname) 
     657      lat_name = 'nav_lat_'//TRIM(phname) 
     658    ENDIF 
     659    lonbound_name = TRIM(lon_name)//'_bounds' 
     660    latbound_name = TRIM(lat_name)//'_bounds' 
    677661  ENDIF 
    678662!- 
     
    682666  W_F(idf)%n_hax = phid 
    683667  W_F(idf)%hax_name(phid,1:2) = (/ lon_name,lat_name /) 
    684   tmp_title = phtitle 
    685668!- 
    686669! 2.0 Longitude 
    687670!- 
    688   IF (l_dbg) WRITE(*,*) "histhori_regular 2.0" 
    689 !- 
    690   IF (rectilinear) THEN 
    691     ndim = 1 
    692     dims(1:1) = (/ W_F(idf)%xid /) 
    693   ENDIF 
    694   iret = NF90_DEF_VAR (nfid,lon_name,NF90_REAL4,dims(1:ndim),nlonid) 
    695   IF (rectilinear) THEN 
    696     iret = NF90_PUT_ATT (nfid,nlonid,'axis',"X") 
    697   ENDIF 
    698   iret = NF90_PUT_ATT (nfid,nlonid,'standard_name',"longitude") 
    699   iret = NF90_PUT_ATT (nfid,nlonid,'units',"degrees_east") 
    700   iret = NF90_PUT_ATT (nfid,nlonid,'valid_min', & 
    701  &                     REAL(MINVAL(plon),KIND=4)) 
    702   iret = NF90_PUT_ATT (nfid,nlonid,'valid_max', & 
    703  &                     REAL(MAXVAL(plon),KIND=4)) 
    704   iret = NF90_PUT_ATT (nfid,nlonid,'long_name',"Longitude") 
    705   iret = NF90_PUT_ATT (nfid,nlonid,'nav_model',TRIM(tmp_title)) 
     671  IF (l_dbg) WRITE(*,*) c_nam//" 2.0" 
     672!- 
     673  i_s = 1; 
     674  IF ( (k_typ == 1).OR.(k_typ == 3) ) THEN 
     675    i_e = 1; wmn = MINVAL(x_1d); wmx = MAXVAL(x_1d); 
     676  ELSEIF (k_typ == 2) THEN 
     677    i_e = 2; wmn = MINVAL(x_2d); wmx = MAXVAL(x_2d); 
     678  ENDIF 
     679  iret = NF90_DEF_VAR(nfid,lon_name,NF90_REAL4,dims(i_s:i_e),nlonid) 
     680  IF (k_typ == 1) THEN 
     681    iret = NF90_PUT_ATT(nfid,nlonid,'axis',"X") 
     682  ENDIF 
     683  iret = NF90_PUT_ATT(nfid,nlonid,'standard_name',"longitude") 
     684  iret = NF90_PUT_ATT(nfid,nlonid,'units',"degrees_east") 
     685  iret = NF90_PUT_ATT(nfid,nlonid,'valid_min',REAL(wmn,KIND=4)) 
     686  iret = NF90_PUT_ATT(nfid,nlonid,'valid_max',REAL(wmx,KIND=4)) 
     687  iret = NF90_PUT_ATT(nfid,nlonid,'long_name',"Longitude") 
     688  iret = NF90_PUT_ATT(nfid,nlonid,'nav_model',TRIM(phtitle)) 
     689!- 
     690  IF (k_typ == 3) THEN 
     691!--- 
     692!-- 2.1 Longitude bounds 
     693!--- 
     694    iret = NF90_PUT_ATT(nfid,nlonid,'bounds',TRIM(lonbound_name)) 
     695    iret = NF90_DEF_VAR(nfid,lonbound_name,NF90_REAL4,dims_b(1:2),nlonidb) 
     696    iret = NF90_PUT_ATT(nfid,nlonidb,'long_name', & 
     697 &          'Boundaries for coordinate variable '//TRIM(lon_name)) 
     698  ENDIF 
    706699!- 
    707700! 3.0 Latitude 
    708701!- 
    709   IF (l_dbg) WRITE(*,*) "histhori_regular 3.0" 
    710 !- 
    711   IF (rectilinear) THEN 
    712     ndim = 1 
    713     dims(1:1) = (/ W_F(idf)%yid /) 
    714   ENDIF 
    715   iret = NF90_DEF_VAR (nfid,lat_name,NF90_REAL4,dims(1:ndim),nlatid) 
    716   IF (rectilinear) THEN 
    717     iret = NF90_PUT_ATT (nfid,nlatid,'axis',"Y") 
    718   ENDIF 
    719   iret = NF90_PUT_ATT (nfid,nlatid,'standard_name',"latitude") 
    720   iret = NF90_PUT_ATT (nfid,nlatid,'units',"degrees_north") 
    721   iret = NF90_PUT_ATT (nfid,nlatid,'valid_min', & 
    722  &                     REAL(MINVAL(plat),KIND=4)) 
    723   iret = NF90_PUT_ATT (nfid,nlatid,'valid_max', & 
    724  &                     REAL(MAXVAL(plat),KIND=4)) 
    725   iret = NF90_PUT_ATT (nfid,nlatid,'long_name',"Latitude") 
    726   iret = NF90_PUT_ATT (nfid,nlatid,'nav_model',TRIM(tmp_title)) 
    727 !- 
    728   iret = NF90_ENDDEF (nfid) 
     702  IF (l_dbg) WRITE(*,*) c_nam//" 3.0" 
     703!- 
     704  i_e = 2; 
     705  IF ( (k_typ == 1).OR.(k_typ == 3) ) THEN 
     706    i_s = 2; wmn = MINVAL(y_1d); wmx = MAXVAL(y_1d); 
     707  ELSEIF (k_typ == 2) THEN 
     708    i_s = 1; wmn = MINVAL(y_2d); wmx = MAXVAL(y_2d); 
     709  ENDIF 
     710  iret = NF90_DEF_VAR(nfid,lat_name,NF90_REAL4,dims(i_s:i_e),nlatid) 
     711  IF (k_typ == 1) THEN 
     712    iret = NF90_PUT_ATT(nfid,nlatid,'axis',"Y") 
     713  ENDIF 
     714!- 
     715  iret = NF90_PUT_ATT(nfid,nlatid,'standard_name',"latitude") 
     716  iret = NF90_PUT_ATT(nfid,nlatid,'units',"degrees_north") 
     717  iret = NF90_PUT_ATT(nfid,nlatid,'valid_min',REAL(wmn,KIND=4)) 
     718  iret = NF90_PUT_ATT(nfid,nlatid,'valid_max',REAL(wmx,KIND=4)) 
     719  iret = NF90_PUT_ATT(nfid,nlatid,'long_name',"Latitude") 
     720  iret = NF90_PUT_ATT(nfid,nlatid,'nav_model',TRIM(phtitle)) 
     721!- 
     722  IF (k_typ == 3) THEN 
     723!--- 
     724!-- 3.1 Latitude bounds 
     725!--- 
     726    iret = NF90_PUT_ATT(nfid,nlatid,'bounds',TRIM(latbound_name)) 
     727    iret = NF90_DEF_VAR(nfid,latbound_name,NF90_REAL4,dims_b(1:2),nlatidb) 
     728    iret = NF90_PUT_ATT(nfid,nlatidb,'long_name', & 
     729 &          'Boundaries for coordinate variable '//TRIM(lat_name)) 
     730  ENDIF 
     731!- 
     732  iret = NF90_ENDDEF(nfid) 
    729733!- 
    730734! 4.0 storing the geographical coordinates 
    731735!- 
    732   IF (l_dbg) WRITE(*,*) "histhori_regular 4.0" 
    733 !- 
    734   orix = W_F(idf)%slab_ori(1) 
    735   oriy = W_F(idf)%slab_ori(2) 
    736   par_szx = W_F(idf)%slab_sz(1) 
    737   par_szy = W_F(idf)%slab_sz(2) 
    738 !- 
    739 ! Transfer the longitude 
    740 !- 
    741   IF (rectilinear) THEN 
    742     iret = NF90_PUT_VAR (nfid,nlonid,plon(orix:orix+par_szx-1,1)) 
    743   ELSE 
    744     iret = NF90_PUT_VAR (nfid,nlonid, & 
    745  &           plon(orix:orix+par_szx-1,oriy:oriy+par_szy-1)) 
    746   ENDIF 
    747 !- 
    748 ! Transfer the latitude 
    749 !- 
    750   IF (rectilinear) THEN 
    751     iret = NF90_PUT_VAR (nfid,nlatid,plat(1,oriy:oriy+par_szy-1)) 
    752   ELSE 
    753     iret = NF90_PUT_VAR (nfid,nlatid, & 
    754  &           plat(orix:orix+par_szx-1,oriy:oriy+par_szy-1)) 
    755   ENDIF 
    756 !- 
    757   iret = NF90_REDEF (nfid) 
    758 !------------------------------ 
    759 END SUBROUTINE histhori_regular 
    760 !=== 
    761 SUBROUTINE histhori_irregular & 
    762  &  (idf,pim,plon,plon_bounds,plat,plat_bounds, & 
    763  &   phname,phtitle,phid) 
    764 !--------------------------------------------------------------------- 
    765 !- This subroutine is made to declare a new horizontale grid. 
    766 !- It has to have the same number of points as 
    767 !- the original and thus in this routine we will only 
    768 !- add two variable (longitude and latitude). 
    769 !- Any variable in the file can thus point to this pair 
    770 !- through an attribute. This routine is very usefull 
    771 !- to allow staggered grids. 
    772 !- 
    773 !- INPUT 
    774 !- 
    775 !- idf         : The id of the file to which the grid should be added 
    776 !- pim         : Size in the longitude direction 
    777 !- plon        : The longitudes 
    778 !- plon_bounds : The boundaries of the grid in longitude 
    779 !- plat        : The latitudes 
    780 !- plat_bounds : Boundaries of the grid in latitude 
    781 !- phname      : The name of grid 
    782 !- phtitle     : The title of the grid 
    783 !- 
    784 !- OUTPUT 
    785 !- 
    786 !- phid    : Id of the created grid 
    787 !--------------------------------------------------------------------- 
    788   IMPLICIT NONE 
    789 !- 
    790   INTEGER,INTENT(IN) :: idf,pim 
    791   REAL,DIMENSION(pim),INTENT(IN) :: plon,plat 
    792   REAL,DIMENSION(:,:),INTENT(IN) :: plon_bounds,plat_bounds 
    793   CHARACTER(LEN=*),INTENT(IN) :: phname,phtitle 
    794   INTEGER,INTENT(OUT) :: phid 
    795 !- 
    796   CHARACTER(LEN=25) :: lon_name,lat_name 
    797   CHARACTER(LEN=30) :: lonbound_name,latbound_name 
    798   CHARACTER(LEN=80) :: tmp_title,tmp_name,longname 
    799   INTEGER :: ndim,dims(2) 
    800   INTEGER :: ndimb,dimsb(2) 
    801   INTEGER :: nbbounds 
    802   INTEGER :: nlonid,nlatid,nlonidb,nlatidb 
    803   INTEGER :: iret,nfid,twoid 
    804   LOGICAL :: transp = .FALSE. 
    805   REAL,ALLOCATABLE,DIMENSION(:,:) :: bounds_trans 
    806   LOGICAL :: l_dbg 
    807 !--------------------------------------------------------------------- 
    808   CALL ipsldbg (old_status=l_dbg) 
    809 !- 
    810 ! 1.0 Check that all fits in the buffers 
    811 !- 
    812   IF (    (pim /= W_F(idf)%full_size(1)) & 
    813  &    .OR.(W_F(idf)%full_size(2) /= 1) ) THEN 
    814     CALL ipslerr (3,"histhori", & 
    815  &   'The new horizontal grid does not have the same size', & 
    816  &   'as the one provided to histbeg. This is not yet ', & 
    817  &   'possible in the hist package.') 
    818   ENDIF 
    819 !- 
    820 ! 1.1 Create all the variables needed 
    821 !- 
    822   IF (l_dbg) WRITE(*,*) 'histhori_irregular 1.0' 
    823 !- 
    824   nfid = W_F(idf)%ncfid 
    825 !- 
    826   IF     (SIZE(plon_bounds,DIM=1) == pim) THEN 
    827     nbbounds = SIZE(plon_bounds,DIM=2) 
    828     transp = .TRUE. 
    829   ELSEIF (SIZE(plon_bounds,DIM=2) == pim) THEN 
    830     nbbounds = SIZE(plon_bounds,DIM=1) 
    831     transp = .FALSE. 
    832   ELSE 
    833     CALL ipslerr (3,"histhori", & 
    834    &  'The boundary variable does not have any axis corresponding', & 
    835    &  'to the size of the longitude or latitude variable', & 
    836    &  '.') 
    837   ENDIF 
    838 !- 
    839   ALLOCATE(bounds_trans(nbbounds,pim)) 
    840 !- 
    841   iret = NF90_DEF_DIM (nfid,'nbnd',nbbounds,twoid) 
    842   ndim = 1 
    843   dims(1) = W_F(idf)%xid 
    844   ndimb = 2 
    845   dimsb(1:2) = (/ twoid,W_F(idf)%xid /) 
    846 !- 
    847   tmp_name = phname 
    848   IF (W_F(idf)%n_hax == 0) THEN 
    849     lon_name = 'nav_lon' 
    850     lat_name = 'nav_lat' 
    851   ELSE 
    852     lon_name = 'nav_lon_'//TRIM(tmp_name) 
    853     lat_name = 'nav_lat_'//TRIM(tmp_name) 
    854   ENDIF 
    855   lonbound_name = TRIM(lon_name)//'_bounds' 
    856   latbound_name = TRIM(lat_name)//'_bounds' 
    857 !- 
    858 ! 1.2 Save the informations 
    859 !- 
    860   phid = W_F(idf)%n_hax+1 
    861   W_F(idf)%n_hax = phid 
    862   W_F(idf)%hax_name(phid,1:2) = (/ lon_name,lat_name /) 
    863   tmp_title = phtitle 
    864 !- 
    865 ! 2.0 Longitude 
    866 !- 
    867   IF (l_dbg) WRITE(*,*) "histhori_irregular 2.0" 
    868 !- 
    869   iret = NF90_DEF_VAR (nfid,lon_name,NF90_REAL4,dims(1:ndim),nlonid) 
    870   iret = NF90_PUT_ATT (nfid,nlonid,'standard_name',"longitude") 
    871   iret = NF90_PUT_ATT (nfid,nlonid,'units',"degrees_east") 
    872   iret = NF90_PUT_ATT (nfid,nlonid,'valid_min', & 
    873  &                     REAL(MINVAL(plon),KIND=4)) 
    874   iret = NF90_PUT_ATT (nfid,nlonid,'valid_max', & 
    875  &                     REAL(MAXVAL(plon),KIND=4)) 
    876   iret = NF90_PUT_ATT (nfid,nlonid,'long_name',"Longitude") 
    877   iret = NF90_PUT_ATT (nfid,nlonid,'nav_model',TRIM(tmp_title)) 
    878 !- 
    879 ! 2.1 Longitude bounds 
    880 !- 
    881   iret = NF90_PUT_ATT (nfid,nlonid,'bounds',TRIM(lonbound_name)) 
    882   iret = NF90_DEF_VAR (nfid,lonbound_name,NF90_REAL4, & 
    883  &                     dimsb(1:ndimb),nlonidb) 
    884   longname = 'Boundaries for coordinate variable '//TRIM(lon_name) 
    885   iret = NF90_PUT_ATT (nfid,nlonidb,'long_name',TRIM(longname)) 
    886 !- 
    887 ! 3.0 Latitude 
    888 !- 
    889   IF (l_dbg) WRITE(*,*) "histhori_irregular 3.0" 
    890 !- 
    891   iret = NF90_DEF_VAR (nfid,lat_name,NF90_REAL4,dims(1:ndim),nlatid) 
    892   iret = NF90_PUT_ATT (nfid,nlatid,'standard_name',"latitude") 
    893   iret = NF90_PUT_ATT (nfid,nlatid,'units',"degrees_north") 
    894   iret = NF90_PUT_ATT (nfid,nlatid,'valid_min', & 
    895  &                     REAL(MINVAL(plat),KIND=4)) 
    896   iret = NF90_PUT_ATT (nfid,nlatid,'valid_max', & 
    897  &                     REAL(MAXVAL(plat),KIND=4)) 
    898   iret = NF90_PUT_ATT (nfid,nlatid,'long_name',"Latitude") 
    899   iret = NF90_PUT_ATT (nfid,nlatid,'nav_model',TRIM(tmp_title)) 
    900 !- 
    901 ! 3.1 Latitude bounds 
    902 !- 
    903   iret = NF90_PUT_ATT (nfid,nlatid,'bounds',TRIM(latbound_name)) 
    904   iret = NF90_DEF_VAR (nfid,latbound_name,NF90_REAL4, & 
    905  &                     dimsb(1:ndimb),nlatidb) 
    906   longname = 'Boundaries for coordinate variable '//TRIM(lat_name) 
    907   iret = NF90_PUT_ATT (nfid,nlatidb,'long_name',TRIM(longname)) 
    908 !- 
    909   iret = NF90_ENDDEF (nfid) 
    910 !- 
    911 ! 4.0 storing the geographical coordinates 
    912 !- 
    913   IF (l_dbg) WRITE(*,*) "histhori_irregular 4.0" 
    914 !- 
    915 ! 4.1 Write the longitude 
    916 !- 
    917   iret = NF90_PUT_VAR (nfid,nlonid,plon(1:pim)) 
    918 !- 
    919 ! 4.2 Write the longitude bounds 
    920 !- 
    921   IF (transp) THEN 
    922     bounds_trans = TRANSPOSE(plon_bounds) 
    923   ELSE 
    924     bounds_trans = plon_bounds 
    925   ENDIF 
    926   iret = NF90_PUT_VAR (nfid,nlonidb,bounds_trans(1:nbbounds,1:pim)) 
    927 !- 
    928 ! 4.3 Write the latitude 
    929 !- 
    930   iret = NF90_PUT_VAR (nfid,nlatid,plat(1:pim)) 
    931 !- 
    932 ! 4.4 Write the latitude bounds 
    933 !- 
    934   IF (transp) THEN 
    935     bounds_trans = TRANSPOSE(plat_bounds) 
    936   ELSE 
    937     bounds_trans = plat_bounds 
    938   ENDIF 
    939   iret = NF90_PUT_VAR (nfid,nlatidb,bounds_trans(1:nbbounds,1:pim)) 
    940 !- 
    941   DEALLOCATE(bounds_trans) 
    942 !- 
    943   iret = NF90_REDEF (nfid) 
    944 !-------------------------------- 
    945 END SUBROUTINE histhori_irregular 
     736  IF (l_dbg) WRITE(*,*) c_nam//" 4.0" 
     737!- 
     738  IF ( (k_typ == 1).OR.(k_typ == 2) ) THEN 
     739    o_x = W_F(idf)%slab_ori(1) 
     740    o_y = W_F(idf)%slab_ori(2) 
     741    s_x = W_F(idf)%slab_siz(1) 
     742    s_y = W_F(idf)%slab_siz(2) 
     743!--- 
     744!-- Transfer the longitude and  the latitude 
     745!--- 
     746    IF     (k_typ == 1) THEN 
     747      iret = NF90_PUT_VAR(nfid,nlonid,x_1d(o_x:o_x+s_x-1)) 
     748      iret = NF90_PUT_VAR(nfid,nlatid,y_1d(o_y:o_y+s_y-1)) 
     749    ELSEIF (k_typ == 2) THEN 
     750      iret = NF90_PUT_VAR(nfid,nlonid, & 
     751 &            x_2d(o_x:o_x+s_x-1,o_y:o_y+s_y-1)) 
     752      iret = NF90_PUT_VAR(nfid,nlatid, & 
     753 &            y_2d(o_x:o_x+s_x-1,o_y:o_y+s_y-1)) 
     754    ENDIF 
     755  ELSEIF (k_typ == 3) THEN 
     756!--- 
     757!-- Transfer the longitude and the longitude bounds 
     758!--- 
     759    iret = NF90_PUT_VAR(nfid,nlonid,x_1d(1:pim)) 
     760!--- 
     761    IF (transp) THEN 
     762      bounds_trans = TRANSPOSE(x_bnds) 
     763    ELSE 
     764      bounds_trans = x_bnds 
     765    ENDIF 
     766    iret = NF90_PUT_VAR(nfid,nlonidb,bounds_trans(1:nbbounds,1:pim)) 
     767!--- 
     768!-- Transfer the latitude and the latitude bounds 
     769!--- 
     770    iret = NF90_PUT_VAR(nfid,nlatid,y_1d(1:pim)) 
     771!--- 
     772    IF (transp) THEN 
     773      bounds_trans = TRANSPOSE(y_bnds) 
     774    ELSE 
     775      bounds_trans = y_bnds 
     776    ENDIF 
     777    iret = NF90_PUT_VAR(nfid,nlatidb,bounds_trans(1:nbbounds,1:pim)) 
     778!--- 
     779    DEALLOCATE(bounds_trans) 
     780  ENDIF 
     781!- 
     782  iret = NF90_REDEF(nfid) 
     783!----------------------- 
     784END SUBROUTINE histh_all 
    946785!=== 
    947786SUBROUTINE histvert (idf,pzaxname,pzaxtitle,pzaxunit, & 
     
    12391078 &  (/ W_F(idf)%slab_ori(1),W_F(idf)%slab_ori(2),par_oriz /) 
    12401079  W_F(idf)%W_V(iv)%zsize(1:3) = & 
    1241  &  (/ W_F(idf)%slab_sz(1),W_F(idf)%slab_sz(2),par_szz /) 
     1080 &  (/ W_F(idf)%slab_siz(1),W_F(idf)%slab_siz(2),par_szz /) 
    12421081!- 
    12431082! Is the size of the full array the same as that of the coordinates ? 
     
    12571096! Is the size of the zoom smaller than the coordinates ? 
    12581097!- 
    1259   IF (    (W_F(idf)%full_size(1) < W_F(idf)%slab_sz(1)) & 
    1260  &    .OR.(W_F(idf)%full_size(2) < W_F(idf)%slab_sz(2)) ) THEN 
     1098  IF (    (W_F(idf)%full_size(1) < W_F(idf)%slab_siz(1)) & 
     1099 &    .OR.(W_F(idf)%full_size(2) < W_F(idf)%slab_siz(2)) ) THEN 
    12611100    str70 = & 
    12621101 &   "Size of variable should be greater or equal to those of the zoom" 
    12631102    WRITE(str71,'("Size of XY zoom :",2I4)') & 
    1264  &   W_F(idf)%slab_sz(1),W_F(idf)%slab_sz(2) 
     1103 &   W_F(idf)%slab_siz(1),W_F(idf)%slab_siz(2) 
    12651104    WRITE(str72,'("Size declared for variable ",A," :",2I4)') & 
    12661105 &   TRIM(tmp_name),pxsize,pysize 
     
    18201659!     Thanks Marine for showing us what errors users can make ! 
    18211660!- 
    1822   IF ( (idf < 1).OR.(idf > nb_files) ) THEN 
     1661  IF ( (idf < 1).OR.(idf > nb_files_max) ) THEN 
    18231662    CALL ipslerr (3,"histwrite", & 
    18241663 &    'Illegal file ID in the histwrite of variable',pvarname,' ') 
     
    24152254  INTEGER :: ifile,iret,i_s,i_e 
    24162255!- 
    2417   LOGICAL :: file_exists 
    24182256  LOGICAL :: l_dbg 
    24192257!--------------------------------------------------------------------- 
     
    24252263!- 
    24262264  IF (PRESENT(file)) THEN 
    2427     IF ( (file >= 1).AND.(file <= nb_files) ) THEN 
     2265    IF ( (file >= 1).AND.(file <= nb_files_max) ) THEN 
    24282266      IF (W_F(ifile)%ncfid > 0) THEN 
    24292267        i_s = file 
     
    24402278  ELSE 
    24412279    i_s = 1 
    2442     i_e = nb_files 
     2280    i_e = nb_files_max 
    24432281  ENDIF 
    24442282!- 
     
    24802318!- 
    24812319  IF (PRESENT(idf)) THEN 
    2482     IF ( (idf >= 1).AND.(idf <= nb_files) ) THEN 
     2320    IF ( (idf >= 1).AND.(idf <= nb_files_max) ) THEN 
    24832321      IF (W_F(ifile)%ncfid > 0) THEN 
    24842322        i_s = idf 
     
    24952333  ELSE 
    24962334    i_s = 1 
    2497     i_e = nb_files 
     2335    i_e = nb_files_max 
    24982336  ENDIF 
    24992337!- 
Note: See TracChangeset for help on using the changeset viewer.