Changeset 962 for IOIPSL/trunk/src/histcom.f90
- Timestamp:
- 03/30/10 11:36:17 (14 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
IOIPSL/trunk/src/histcom.f90
r957 r962 37 37 !- 38 38 INTERFACE histbeg 39 MODULE PROCEDURE histb eg_totreg,histbeg_regular,histbeg_irregular39 MODULE PROCEDURE histb_reg1d,histb_reg2d,histb_irreg 40 40 END INTERFACE 41 41 !- 42 42 INTERFACE histhori 43 MODULE PROCEDURE histh ori_regular,histhori_irregular43 MODULE PROCEDURE histh_reg1d,histh_reg2d,histh_irreg 44 44 END INTERFACE 45 45 !- … … 129 129 INTEGER :: tid,xid,yid 130 130 !-General definitions in the NETCDF file 131 INTEGER,DIMENSION(2) :: full_size=0,slab_ori,slab_s z131 INTEGER,DIMENSION(2) :: full_size=0,slab_ori,slab_siz 132 132 !-The horizontal axes 133 133 CHARACTER(LEN=25),DIMENSION(nb_hax_max,2) :: hax_name … … 145 145 TYPE(T_D_F),DIMENSION(nb_files_max),SAVE :: W_F 146 146 !- 147 ! Counter of elements148 !-149 INTEGER,SAVE :: nb_files=0150 !-151 147 ! A list of functions which require special action 152 148 ! (Needs to be updated when functions are added … … 162 158 !=== 163 159 !- 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 !- 160 SUBROUTINE 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) 215 166 !--------------------------------------------------------------------- 216 167 IMPLICIT NONE … … 224 175 REAL,INTENT(IN) :: pdate0,pdeltat 225 176 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 !------------------------- 186 END SUBROUTINE histb_reg1d 249 187 !=== 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 !- 188 SUBROUTINE 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) 304 194 !--------------------------------------------------------------------- 305 195 IMPLICIT NONE … … 312 202 REAL,INTENT(IN) :: pdate0,pdeltat 313 203 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 !------------------------- 213 END SUBROUTINE histb_reg2d 428 214 !=== 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 !- 215 SUBROUTINE 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) 472 220 !--------------------------------------------------------------------- 473 221 IMPLICIT NONE … … 480 228 REAL,INTENT(IN) :: pdate0,pdeltat 481 229 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 !------------------------- 238 END SUBROUTINE histb_irreg 239 !=== 240 SUBROUTINE 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 483 320 !- 484 321 INTEGER :: nfid,iret,m_c 485 322 CHARACTER(LEN=120) :: file 486 323 CHARACTER(LEN=30) :: timenow 324 CHARACTER(LEN=11) :: c_nam 487 325 LOGICAL :: l_dbg 488 326 !--------------------------------------------------------------------- 489 327 CALL ipsldbg (old_status=l_dbg) 490 328 !- 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 495 336 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 500 355 !- 501 356 ! 1.0 Transfering into the common for future use 502 357 !- 503 IF (l_dbg) WRITE(*,*) "histbeg_irregular1.0"358 IF (l_dbg) WRITE(*,*) c_nam//" 1.0" 504 359 !- 505 360 W_F(idf)%itau0 = pitau0 … … 509 364 ! 2.0 Initializes all variables for this file 510 365 !- 511 IF (l_dbg) WRITE(*,*) "histbeg_irregular2.0"366 IF (l_dbg) WRITE(*,*) c_nam//" 2.0" 512 367 !- 513 368 W_F(idf)%n_var = 0 … … 516 371 W_F(idf)%n_zax = 0 517 372 !- 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 520 380 !- 521 381 ! 3.0 Opening netcdf file and defining dimensions 522 382 !- 523 IF (l_dbg) WRITE(*,*) "histbeg_irregular3.0"383 IF (l_dbg) WRITE(*,*) c_nam//" 3.0" 524 384 !- 525 385 ! Add DOMAIN number and ".nc" suffix in file name if needed 526 386 !- 527 file = pfilename387 file = nc_name 528 388 CALL flio_dom_file (file,domain_id) 529 389 !- 530 390 ! 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') 534 396 m_c = NF90_CLOBBER 535 CASE( 64)397 CASE('64') 536 398 m_c = IOR(NF90_CLOBBER,NF90_64BIT_OFFSET) 537 399 CASE DEFAULT 538 400 CALL ipslerr (3,"histbeg", & 539 & 'Invalid argument nb_bitsfor file :',TRIM(file), &401 & 'Invalid argument mode for file :',TRIM(file), & 540 402 & 'Supported values are 32 or 64') 541 403 END SELECT … … 545 407 !- 546 408 ! 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 551 422 !- 552 423 ! 4.0 Declaring the geographical coordinates and other attributes 553 424 !- 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)) 561 430 lock_modname = .TRUE. 562 431 CALL ioget_timestamp (timenow) 563 iret = NF90_PUT_ATT 432 iret = NF90_PUT_ATT(nfid,NF90_GLOBAL,'TimeStamp',TRIM(timenow)) 564 433 !- 565 434 ! 5.0 Saving some important information on this file in the common 566 435 !- 567 IF (l_dbg) WRITE(*,*) "histbeg_irregular5.0"436 IF (l_dbg) WRITE(*,*) c_nam//" 5.0" 568 437 !- 569 438 IF (PRESENT(domain_id)) THEN … … 571 440 ENDIF 572 441 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 574 449 !- 575 450 ! 6.0 storing the geographical coordinates 576 451 !- 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 !----------------------- 466 END SUBROUTINE histb_all 584 467 !=== 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. 468 SUBROUTINE 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 !------------------------- 483 END SUBROUTINE histh_reg1d 484 !=== 485 SUBROUTINE 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 !------------------------- 500 END SUBROUTINE histh_reg2d 501 !=== 502 SUBROUTINE 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 !------------------------- 519 END SUBROUTINE histh_irreg 520 !=== 521 SUBROUTINE 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. 589 527 !- It has to have the same number of points as 590 528 !- the original and thus in this routine we will only … … 596 534 !- INPUT 597 535 !- 536 !- k_typ : Type of the grid (1 rectilinear, 2 regular, 3 irregular) 598 537 !- idf : The id of the file to which the grid should be added 599 538 !- 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) 603 540 !- phname : The name of grid 604 541 !- phtitle : The title of the grid … … 608 545 !- phid : Id of the created grid 609 546 !- 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 614 560 !--------------------------------------------------------------------- 615 561 IMPLICIT NONE 616 562 !- 563 INTEGER,INTENT(IN) :: k_typ 617 564 INTEGER,INTENT(IN) :: idf,pim,pjm 618 REAL,INTENT(IN),DIMENSION(pim,pjm) :: plon,plat619 565 CHARACTER(LEN=*),INTENT(IN) :: phname,phtitle 620 566 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 622 570 !- 623 571 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 627 580 INTEGER :: nlonid,nlatid 628 INTEGER :: o rix,oriy,par_szx,par_szy581 INTEGER :: o_x,o_y,s_x,s_y 629 582 INTEGER :: iret,nfid 630 LOGICAL :: rectilinear583 CHARACTER(LEN=11) :: c_nam 631 584 LOGICAL :: l_dbg 632 585 !--------------------------------------------------------------------- 633 586 CALL ipsldbg (old_status=l_dbg) 634 587 !- 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 !- 635 599 ! 1.0 Check that all fits in the buffers 636 600 !- 637 601 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 639 604 CALL ipslerr (3,"histhori", & 640 605 & 'The new horizontal grid does not have the same size', & … … 643 608 ENDIF 644 609 !- 645 IF (PRESENT(opt_rectilinear)) THEN646 rectilinear = opt_rectilinear647 ELSE648 rectilinear = .FALSE.649 ENDIF650 !-651 610 ! 1.1 Create all the variables needed 652 611 !- 653 IF (l_dbg) WRITE(*,*) "histhori_regular1.0"612 IF (l_dbg) WRITE(*,*) c_nam//" 1.0" 654 613 !- 655 614 nfid = W_F(idf)%ncfid 656 615 !- 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 !- 658 633 dims(1:2) = (/ W_F(idf)%xid,W_F(idf)%yid /) 659 634 !- 660 tmp_name = phname 661 IF (rectilinear) THEN 635 IF (k_typ == 1) THEN 662 636 IF (W_F(idf)%n_hax == 0) THEN 663 637 lon_name = 'lon' 664 638 lat_name = 'lat' 665 639 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 670 644 IF (W_F(idf)%n_hax == 0) THEN 671 645 lon_name = 'nav_lon' 672 646 lat_name = 'nav_lat' 673 647 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' 677 661 ENDIF 678 662 !- … … 682 666 W_F(idf)%n_hax = phid 683 667 W_F(idf)%hax_name(phid,1:2) = (/ lon_name,lat_name /) 684 tmp_title = phtitle685 668 !- 686 669 ! 2.0 Longitude 687 670 !- 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 706 699 !- 707 700 ! 3.0 Latitude 708 701 !- 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) 729 733 !- 730 734 ! 4.0 storing the geographical coordinates 731 735 !- 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 !----------------------- 784 END SUBROUTINE histh_all 946 785 !=== 947 786 SUBROUTINE histvert (idf,pzaxname,pzaxtitle,pzaxunit, & … … 1239 1078 & (/ W_F(idf)%slab_ori(1),W_F(idf)%slab_ori(2),par_oriz /) 1240 1079 W_F(idf)%W_V(iv)%zsize(1:3) = & 1241 & (/ W_F(idf)%slab_s z(1),W_F(idf)%slab_sz(2),par_szz /)1080 & (/ W_F(idf)%slab_siz(1),W_F(idf)%slab_siz(2),par_szz /) 1242 1081 !- 1243 1082 ! Is the size of the full array the same as that of the coordinates ? … … 1257 1096 ! Is the size of the zoom smaller than the coordinates ? 1258 1097 !- 1259 IF ( (W_F(idf)%full_size(1) < W_F(idf)%slab_s z(1)) &1260 & .OR.(W_F(idf)%full_size(2) < W_F(idf)%slab_s z(2)) ) THEN1098 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 1261 1100 str70 = & 1262 1101 & "Size of variable should be greater or equal to those of the zoom" 1263 1102 WRITE(str71,'("Size of XY zoom :",2I4)') & 1264 & W_F(idf)%slab_s z(1),W_F(idf)%slab_sz(2)1103 & W_F(idf)%slab_siz(1),W_F(idf)%slab_siz(2) 1265 1104 WRITE(str72,'("Size declared for variable ",A," :",2I4)') & 1266 1105 & TRIM(tmp_name),pxsize,pysize … … 1820 1659 ! Thanks Marine for showing us what errors users can make ! 1821 1660 !- 1822 IF ( (idf < 1).OR.(idf > nb_files ) ) THEN1661 IF ( (idf < 1).OR.(idf > nb_files_max) ) THEN 1823 1662 CALL ipslerr (3,"histwrite", & 1824 1663 & 'Illegal file ID in the histwrite of variable',pvarname,' ') … … 2415 2254 INTEGER :: ifile,iret,i_s,i_e 2416 2255 !- 2417 LOGICAL :: file_exists2418 2256 LOGICAL :: l_dbg 2419 2257 !--------------------------------------------------------------------- … … 2425 2263 !- 2426 2264 IF (PRESENT(file)) THEN 2427 IF ( (file >= 1).AND.(file <= nb_files ) ) THEN2265 IF ( (file >= 1).AND.(file <= nb_files_max) ) THEN 2428 2266 IF (W_F(ifile)%ncfid > 0) THEN 2429 2267 i_s = file … … 2440 2278 ELSE 2441 2279 i_s = 1 2442 i_e = nb_files 2280 i_e = nb_files_max 2443 2281 ENDIF 2444 2282 !- … … 2480 2318 !- 2481 2319 IF (PRESENT(idf)) THEN 2482 IF ( (idf >= 1).AND.(idf <= nb_files ) ) THEN2320 IF ( (idf >= 1).AND.(idf <= nb_files_max) ) THEN 2483 2321 IF (W_F(ifile)%ncfid > 0) THEN 2484 2322 i_s = idf … … 2495 2333 ELSE 2496 2334 i_s = 1 2497 i_e = nb_files 2335 i_e = nb_files_max 2498 2336 ENDIF 2499 2337 !-
Note: See TracChangeset
for help on using the changeset viewer.