Changeset 5600 for branches/2014/dev_r4650_UKMO14.12_STAND_ALONE_OBSOPER/NEMOGCM/TOOLS/SIREN/src/merge_bathy.f90
- Timestamp:
- 2015-07-15T17:46:12+02:00 (9 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2014/dev_r4650_UKMO14.12_STAND_ALONE_OBSOPER/NEMOGCM/TOOLS/SIREN/src/merge_bathy.f90
r4213 r5600 7 7 ! 8 8 ! DESCRIPTION: 9 !> @file 9 10 !> @brief 10 11 !> This program merge bathymetry file at boundaries. 11 12 !> 12 13 !> @details 14 !> @section sec1 method 13 15 !> Coarse grid Bathymetry is interpolated on fine grid. 14 !> Then fine Bathymetry and refined coarse bathymetry are merged at boundaries. 15 !> 16 !> BathyFine= weight * BathyCoarse + (1-weight)*BathyFine 17 !> 18 !> The weight function used is : 0.5 + 0.5*COS( (pi*dist) / width ) 19 !> 20 !> @author 21 !> J.Paul 16 !> Then fine Bathymetry and refined coarse bathymetry are merged at boundaries.<br/> 17 !> @f[BathyFine= Weight * BathyCoarse + (1-Weight)*BathyFine@f] 18 !> The weight function used is :<br/> 19 !> @f[Weight = 0.5 + 0.5*COS( \frac{\pi*dist}{width} )@f]<br/> 20 !> with 21 !> - dist : number of point to border 22 !> - width : boundary size 23 !> 24 !> @section sec2 how to 25 !> to merge bathymetry file:<br/> 26 !> @code{.sh} 27 !> ./SIREN/bin/merge_bathy merge_bathy.nam 28 !> @endcode 29 !> 30 !> merge_bathy.nam comprise 8 namelists: 31 !> - logger namelist (namlog) 32 !> - config namelist (namcfg) 33 !> - coarse grid namelist (namcrs) 34 !> - fine grid namelist (namfin) 35 !> - variable namelist (namvar) 36 !> - nesting namelist (namnst) 37 !> - boundary namelist (nambdy) 38 !> - output namelist (namout) 39 !> 40 !> @note 41 !> All namelists have to be in file merge_bathy.nam, 42 !> however variables of those namelists are all optional. 43 !> 44 !> * _logger namelist (namlog)_: 45 !> - cn_logfile : logger filename 46 !> - cn_verbosity : verbosity ('trace','debug','info', 47 !> 'warning','error','fatal') 48 !> - in_maxerror : maximum number of error allowed 49 !> 50 !> * _config namelist (namcfg)_: 51 !> - cn_varcfg : variable configuration file (see ./SIREN/cfg/variable.cfg) 52 !> 53 !> * _coarse grid namelist (namcrs)_: 54 !> - cn_bathy0 : bathymetry file 55 !> - in_perio0 : NEMO periodicity index (see Model Boundary Condition in 56 !> [NEMO documentation](http://www.nemo-ocean.eu/About-NEMO/Reference-manuals)) 57 !> 58 !> * _fine grid namelist (namfin)_: 59 !> - cn_bathy1 : bathymetry file 60 !> - in_perio1 : NEMO periodicity index 61 !> 62 !> * _variable namelist (namvar)_: 63 !> - cn_varinfo : list of variable and extra information about request(s) 64 !> to be used.<br/> 65 !> each elements of *cn_varinfo* is a string character.<br/> 66 !> it is composed of the variable name follow by ':', 67 !> then request(s) to be used on this variable.<br/> 68 !> request could be: 69 !> - interpolation method 70 !> 71 !> requests must be separated by ';'.<br/> 72 !> order of requests does not matter.<br/> 73 !> 74 !> informations about available method could be find in 75 !> @ref interp modules.<br/> 76 !> Example: 'bathymetry: cubic' 77 !> @note 78 !> If you do not specify a method which is required, 79 !> default one is apply. 80 !> @warning 81 !> variable name must be __Bathymetry__ here. 82 !> 83 !> * _nesting namelist (namnst)_: 84 !> - in_rhoi : refinement factor in i-direction 85 !> - in_rhoj : refinement factor in j-direction 86 !> 87 !> * _boundary namelist (nambdy)_: 88 !> - ln_north : use north boundary or not 89 !> - ln_south : use south boundary or not 90 !> - ln_east : use east boundary or not 91 !> - ln_west : use west boundary or not 92 !> - cn_north : north boundary indices on fine grid<br/> 93 !> *cn_north* is a string character defining boundary 94 !> segmentation.<br/> 95 !> segments are separated by '|'.<br/> 96 !> each segments of the boundary is composed of: 97 !> - orthogonal indice (.ie. for north boundary, 98 !> J-indice where boundary are). 99 !> - first indice of boundary (I-indice for north boundary) 100 !> - last indice of boundary (I-indice for north boundary)<br/> 101 !> indices must be separated by ',' .<br/> 102 !> - optionally, boundary size could be added between '(' and ')' 103 !> in the first segment defined. 104 !> @note 105 !> boundary size is the same for all segments of one boundary. 106 !> 107 !> Examples: 108 !> - cn_north='index1,first1,last1(width)' 109 !> - cn_north='index1(width),first1,last1|index2,first2,last2' 110 !> - cn_south : south boundary indices on fine grid<br/> 111 !> - cn_east : east boundary indices on fine grid<br/> 112 !> - cn_west : west boundary indices on fine grid<br/> 113 !> - ln_oneseg: use only one segment for each boundary or not 114 !> 115 !> * _output namelist (namout)_: 116 !> - cn_fileout : merged bathymetry file 117 !> 118 !> @author J.Paul 22 119 ! REVISION HISTORY: 23 !> @date Nov, 2013 - Initial Version 24 ! 120 !> @date November, 2013 - Initial Version 121 !> @date Sepember, 2014 122 !> - add header for user 123 !> 25 124 !> @note Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 26 !>27 !> @todo28 125 !---------------------------------------------------------------------- 29 !> @code30 126 PROGRAM merge_bathy 31 127 … … 43 139 USE boundary ! boundary manager 44 140 USE iom ! I/O manager 45 USE dom ! domain manager46 141 USE grid ! grid manager 47 142 USE extrap ! extrapolation manager … … 49 144 USE filter ! filter manager 50 145 USE mpp ! MPP manager 146 USE dom ! domain manager 51 147 USE iom_mpp ! MPP I/O manager 148 USE iom_dom ! DOM I/O manager 52 149 53 150 IMPLICIT NONE … … 60 157 INTEGER(i4) :: il_status 61 158 INTEGER(i4) :: il_fileid 62 INTEGER(i4) :: il_atti d159 INTEGER(i4) :: il_attind 63 160 INTEGER(i4) :: il_imin0 64 161 INTEGER(i4) :: il_imax0 … … 66 163 INTEGER(i4) :: il_jmax0 67 164 INTEGER(i4) , DIMENSION(ip_maxdim) :: il_rho 68 INTEGER(i4) , DIMENSION(2,2) :: il_offset 69 INTEGER(i4) , DIMENSION(2,2,2) :: il_ind 70 ! INTEGER(i4) , DIMENSION(:,:,:,:), ALLOCATABLE :: il_value 165 INTEGER(i4) , DIMENSION(2,2) :: il_ind 71 166 72 167 LOGICAL :: ll_exist … … 76 171 REAL(dp) , DIMENSION(:,:,:,:), ALLOCATABLE :: dl_weight 77 172 78 TYPE(T FILE):: tl_bathy079 TYPE(T FILE):: tl_bathy1173 TYPE(TMPP) :: tl_bathy0 174 TYPE(TMPP) :: tl_bathy1 80 175 TYPE(TFILE) :: tl_fileout 81 176 … … 85 180 TYPE(TVAR) :: tl_lon 86 181 TYPE(TVAR) :: tl_lat 87 ! TYPE(TVAR) :: tl_depth88 ! TYPE(TVAR) :: tl_time89 182 90 183 TYPE(TDIM) , DIMENSION(ip_maxdim) :: tl_dim … … 98 191 99 192 ! namelist variable 193 ! namlog 100 194 CHARACTER(LEN=lc) :: cn_logfile = 'merge_bathy.log' 101 195 CHARACTER(LEN=lc) :: cn_verbosity = 'warning' 102 196 INTEGER(i4) :: in_maxerror = 5 197 198 ! namcfg 199 CHARACTER(LEN=lc) :: cn_varcfg = 'variable.cfg' 200 201 ! namcrs 103 202 CHARACTER(LEN=lc) :: cn_bathy0 = '' 104 203 INTEGER(i4) :: in_perio0 = -1 105 204 205 ! namfin 106 206 CHARACTER(LEN=lc) :: cn_bathy1 = '' 107 207 INTEGER(i4) :: in_perio1 = -1 108 208 109 CHARACTER(LEN=lc) :: cn_varcfg = 'variable.cfg' 110 111 CHARACTER(LEN=lc), DIMENSION(ig_maxvar) :: cn_varinfo = '' 112 113 INTEGER(i4) :: in_imin0 = 0 114 INTEGER(i4) :: in_imax0 = 0 115 INTEGER(i4) :: in_jmin0 = 0 116 INTEGER(i4) :: in_jmax0 = 0 209 ! namvar 210 CHARACTER(LEN=lc), DIMENSION(ip_maxvar) :: cn_varinfo = '' 211 212 ! namnst 117 213 INTEGER(i4) :: in_rhoi = 0 118 214 INTEGER(i4) :: in_rhoj = 0 119 215 216 ! nambdy 120 217 LOGICAL :: ln_north = .TRUE. 121 218 LOGICAL :: ln_south = .TRUE. … … 128 225 CHARACTER(LEN=lc) :: cn_west = '' 129 226 227 ! namout 130 228 CHARACTER(LEN=lc) :: cn_fileout = 'bathy_merged.nc' 131 229 !------------------------------------------------------------------- … … 152 250 153 251 NAMELIST /namnst/ & !< nesting namelist 154 & in_imin0, & !< i-direction lower left point indice on coarse grid155 & in_imax0, & !< i-direction upper right point indice on coarse grid156 & in_jmin0, & !< j-direction lower left point indice on coarse grid157 & in_jmax0, & !< j-direction upper right point indice on coarse grid158 252 & in_rhoi, & !< refinement factor in i-direction 159 253 & in_rhoj !< refinement factor in j-direction … … 171 265 172 266 NAMELIST /namout/ & !< output namelist 173 & cn_fileout !< fine grid merged bathymetry file267 & cn_fileout !< fine grid merged bathymetry file 174 268 !------------------------------------------------------------------- 175 269 176 ! 1-namelist177 ! 1-1get namelist270 ! namelist 271 ! get namelist 178 272 il_narg=COMMAND_ARGUMENT_COUNT() !f03 intrinsec 179 273 IF( il_narg/=1 )THEN … … 184 278 ENDIF 185 279 186 ! 1-2read namelist280 ! read namelist 187 281 INQUIRE(FILE=TRIM(cl_namelist), EXIST=ll_exist) 188 282 IF( ll_exist )THEN … … 203 297 204 298 READ( il_fileid, NML = namlog ) 205 ! 1-2-1define log file206 CALL logger_open(TRIM(cn_logfile),TRIM(cn_verbosity) )299 ! define log file 300 CALL logger_open(TRIM(cn_logfile),TRIM(cn_verbosity), in_maxerror) 207 301 CALL logger_header() 208 302 209 303 READ( il_fileid, NML = namcfg ) 210 ! 1-2-2get variable extra information304 ! get variable extra information 211 305 CALL var_def_extra(TRIM(cn_varcfg)) 212 306 … … 214 308 READ( il_fileid, NML = namfin ) 215 309 READ( il_fileid, NML = namvar ) 216 ! 1-2-3add user change in extra information310 ! add user change in extra information 217 311 CALL var_chg_extra(cn_varinfo) 218 312 … … 234 328 ENDIF 235 329 236 ! 2-open files330 ! open files 237 331 IF( TRIM(cn_bathy0) /= '' )THEN 238 tl_bathy0= file_init(TRIM(cn_bathy0),id_perio=in_perio0)239 CALL iom_open(tl_bathy0)332 tl_bathy0=mpp_init( file_init(TRIM(cn_bathy0)), id_perio=in_perio0) 333 CALL grid_get_info(tl_bathy0) 240 334 ELSE 241 335 CALL logger_fatal("MERGE BATHY: can not find coarse grid bathymetry "//& … … 244 338 245 339 IF( TRIM(cn_bathy1) /= '' )THEN 246 tl_bathy1= file_init(TRIM(cn_bathy1),id_perio=in_perio1)247 CALL iom_open(tl_bathy1)340 tl_bathy1=mpp_init( file_init(TRIM(cn_bathy1)), id_perio=in_perio1) 341 CALL grid_get_info(tl_bathy1) 248 342 ELSE 249 343 CALL logger_fatal("MERGE BATHY: can not find fine grid bathymetry "//& … … 251 345 ENDIF 252 346 253 ! 3-check254 ! 3-1check output file do not already exist347 ! check 348 ! check output file do not already exist 255 349 INQUIRE(FILE=TRIM(cn_fileout), EXIST=ll_exist) 256 350 IF( ll_exist )THEN … … 259 353 ENDIF 260 354 261 ! 3-2check namelist262 ! 3-2-1check refinament factor355 ! check namelist 356 ! check refinament factor 263 357 il_rho(:)=1 264 358 IF( in_rhoi < 1 .OR. in_rhoj < 1 )THEN … … 270 364 ENDIF 271 365 272 !3-2-2 check domain indices 273 IF( in_imin0 < 1 .OR. in_imax0 < 1 .OR. in_jmin0 < 1 .OR. in_jmax0 < 1)THEN 274 ! compute coarse grid indices around fine grid 275 il_ind(:,:,:)=grid_get_coarse_index(tl_bathy0, tl_bathy1 ) 276 277 il_imin0=il_ind(1,1,1) ; il_imax0=il_ind(1,2,1) 278 il_jmin0=il_ind(2,1,1) ; il_jmax0=il_ind(2,2,1) 279 280 il_offset(:,:)=il_ind(:,:,2) 281 ELSE 282 il_imin0=in_imin0 ; il_imax0=in_imax0 283 il_jmin0=in_jmin0 ; il_jmax0=in_jmax0 284 285 il_offset(1,:)=FLOOR(REAL(il_rho(jp_I)-1,dp)*0.5) 286 il_offset(2,:)=FLOOR(REAL(il_rho(jp_J)-1,dp)*0.5) 287 ENDIF 288 289 !3-3 check domain validity 366 ! check domain indices 367 ! compute coarse grid indices around fine grid 368 il_ind(:,:)=grid_get_coarse_index(tl_bathy0, tl_bathy1, & 369 & id_rho=il_rho(:) ) 370 371 il_imin0=il_ind(1,1) ; il_imax0=il_ind(1,2) 372 il_jmin0=il_ind(2,1) ; il_jmax0=il_ind(2,2) 373 374 ! check domain validity 290 375 CALL grid_check_dom(tl_bathy0, il_imin0, il_imax0, il_jmin0, il_jmax0) 291 376 292 ! 3-4check coincidence between coarse and fine grid377 ! check coincidence between coarse and fine grid 293 378 CALL grid_check_coincidence( tl_bathy0, tl_bathy1, & 294 379 & il_imin0, il_imax0, & … … 296 381 & il_rho(:) ) 297 382 298 !4- read or compute boundary 299 tl_var=iom_read_var(tl_bathy1,'Bathymetry') 383 ! open mpp files 384 CALL iom_mpp_open(tl_bathy1) 385 386 ! read or compute boundary 387 tl_var=iom_mpp_read_var(tl_bathy1,'Bathymetry') 388 389 ! close mpp files 390 CALL iom_mpp_close(tl_bathy1) 300 391 301 392 tl_bdy(:)=boundary_init(tl_var, ln_north, ln_south, ln_east, ln_west, & … … 303 394 & ln_oneseg ) 304 395 305 ! 5-get boundary on coarse grid306 ! 5-1 define refined bathymetry table(for coarse grid)396 ! get boundary on coarse grid 397 ! define refined bathymetry array (for coarse grid) 307 398 dl_fill=tl_var%d_fill 308 399 ALLOCATE( dl_refined(tl_var%t_dim(1)%i_len, & … … 313 404 dl_refined(:,:,:,:)=dl_fill 314 405 315 ! 5-2 define weight table406 ! define weight array 316 407 ALLOCATE( dl_weight(tl_var%t_dim(1)%i_len, & 317 408 & tl_var%t_dim(2)%i_len, & … … 320 411 dl_weight(:,:,:,:)=dl_fill 321 412 322 ! 5-3compute coarse grid refined bathymetry on boundary.413 ! compute coarse grid refined bathymetry on boundary. 323 414 DO jk=1,ip_ncard 324 325 415 CALL merge_bathy_get_boundary(tl_bathy0, tl_bathy1, tl_bdy(jk), & 326 416 & il_rho(:), & … … 330 420 ENDDO 331 421 332 ! 6-merge bathy on boundary422 ! merge bathy on boundary 333 423 DO jl=1,tl_var%t_dim(4)%i_len 334 424 DO jk=1,tl_var%t_dim(3)%i_len … … 348 438 DEALLOCATE(dl_refined) 349 439 350 ! 7-create file440 ! create file 351 441 tl_fileout=file_init(TRIM(cn_fileout),id_perio=in_perio1) 352 442 353 ! 7-1add dimension354 tl_dim(:)= tl_var%t_dim(:)443 ! add dimension 444 tl_dim(:)=dim_copy(tl_var%t_dim(:)) 355 445 356 446 DO ji=1,ip_maxdim … … 358 448 ENDDO 359 449 360 ! 7-2add variables450 ! add variables 361 451 IF( ALL( tl_dim(1:2)%l_use ) )THEN 362 363 tl_lon=iom_read_var(tl_bathy1,'longitude') 452 ! open mpp files 453 CALL iom_mpp_open(tl_bathy1) 454 455 ! add longitude 456 tl_lon=iom_mpp_read_var(tl_bathy1,'longitude') 364 457 CALL file_add_var(tl_fileout, tl_lon) 365 458 CALL var_clean(tl_lon) 366 459 367 tl_lat=iom_read_var(tl_bathy1,'latitude') 460 ! add latitude 461 tl_lat=iom_mpp_read_var(tl_bathy1,'latitude') 368 462 CALL file_add_var(tl_fileout, tl_lat) 369 463 CALL var_clean(tl_lat) 370 464 465 ! close mpp files 466 CALL iom_mpp_close(tl_bathy1) 371 467 ENDIF 372 468 … … 375 471 376 472 ! only 2 first dimension to be used 377 tl_dim(:)= tl_fileout%t_dim(:)473 tl_dim(:)=dim_copy(tl_fileout%t_dim(:)) 378 474 tl_dim(3:4)%l_use=.FALSE. 379 475 tl_var=var_init('weight',dl_weight(:,:,:,:),td_dim=tl_dim(:),dd_fill=dl_fill) … … 381 477 CALL var_clean(tl_var) 382 478 383 ! 7-3add some attribute479 ! add some attribute 384 480 tl_att=att_init("Created_by","SIREN merge_bathy") 385 481 CALL file_add_att(tl_fileout, tl_att) … … 395 491 CALL file_add_att(tl_fileout, tl_att) 396 492 397 ! a voir398 493 ! add attribute periodicity 399 il_atti d=0494 il_attind=0 400 495 IF( ASSOCIATED(tl_fileout%t_att) )THEN 401 il_atti d=att_get_id(tl_fileout%t_att(:),'periodicity')402 ENDIF 403 IF( tl_bathy1%i_perio >= 0 .AND. il_atti d == 0 )THEN496 il_attind=att_get_index(tl_fileout%t_att(:),'periodicity') 497 ENDIF 498 IF( tl_bathy1%i_perio >= 0 .AND. il_attind == 0 )THEN 404 499 tl_att=att_init('periodicity',tl_bathy1%i_perio) 405 500 CALL file_add_att(tl_fileout,tl_att) 406 501 ENDIF 407 502 408 il_atti d=0503 il_attind=0 409 504 IF( ASSOCIATED(tl_fileout%t_att) )THEN 410 il_atti d=att_get_id(tl_fileout%t_att(:),'ew_overlap')411 ENDIF 412 IF( tl_bathy1%i_ew >= 0 .AND. il_atti d == 0 )THEN505 il_attind=att_get_index(tl_fileout%t_att(:),'ew_overlap') 506 ENDIF 507 IF( tl_bathy1%i_ew >= 0 .AND. il_attind == 0 )THEN 413 508 tl_att=att_init('ew_overlap',tl_bathy1%i_ew) 414 509 CALL file_add_att(tl_fileout,tl_att) 415 510 ENDIF 416 511 417 ! 7-4create file512 ! create file 418 513 CALL iom_create(tl_fileout) 419 514 420 ! 7-5write file515 ! write file 421 516 CALL iom_write_file(tl_fileout) 422 517 423 ! 7-6close file518 ! close file 424 519 CALL iom_close(tl_fileout) 425 520 426 CALL iom_close(tl_bathy1) 427 CALL iom_close(tl_bathy0) 428 429 !8- clean 521 ! clean 522 CALL att_clean(tl_att) 430 523 CALL file_clean(tl_fileout) 431 CALL file_clean(tl_bathy1)432 CALL file_clean(tl_bathy0)524 CALL mpp_clean(tl_bathy1) 525 CALL mpp_clean(tl_bathy0) 433 526 DEALLOCATE(dl_weight) 434 527 … … 437 530 CALL logger_close() 438 531 439 !> @endcode440 532 CONTAINS 441 533 !------------------------------------------------------------------- 442 534 !> @brief 443 !> This subroutine 535 !> This subroutine compute refined bathymetry on boundary from coarse grid. 444 536 !> 445 !> @details 537 !> @author J.Paul 538 !> @date November, 2013 - Initial Version 446 539 !> 447 !> @author J.Paul 448 !> - Nov, 2013- Initial Version 449 !> 450 !> @param[in] 451 !> @todo 540 !> @param[in] td_bathy0 coarse grid bathymetry file structure 541 !> @param[in] td_bathy1 fine grid bathymetry file structure 542 !> @param[in] td_bdy boundary structure 543 !> @param[in] id_rho array of refinement factor 544 !> @param[inout] dd_refined array of refined bathymetry 545 !> @param[inout] dd_weight array of weight 546 !> @param[in] dd_fill fillValue 452 547 !------------------------------------------------------------------- 453 !> @code454 548 SUBROUTINE merge_bathy_get_boundary( td_bathy0, td_bathy1, td_bdy, & 455 549 & id_rho, & … … 459 553 460 554 ! Argument 461 TYPE(T FILE), INTENT(IN ) :: td_bathy0462 TYPE(T FILE), INTENT(IN ) :: td_bathy1555 TYPE(TMPP) , INTENT(IN ) :: td_bathy0 556 TYPE(TMPP) , INTENT(IN ) :: td_bathy1 463 557 TYPE(TBDY) , INTENT(IN ) :: td_bdy 464 558 INTEGER(i4), DIMENSION(:) , INTENT(IN ) :: id_rho … … 478 572 INTEGER(i4) :: il_jmax0 479 573 480 INTEGER(i4) :: il_imin481 INTEGER(i4) :: il_imax482 INTEGER(i4) :: il_jmin483 INTEGER(i4) :: il_jmax484 485 574 INTEGER(i4), DIMENSION(2,2) :: il_offset 486 INTEGER(i4), DIMENSION(2,2 ,2):: il_ind575 INTEGER(i4), DIMENSION(2,2) :: il_ind 487 576 488 577 REAL(dp) , DIMENSION(:) , ALLOCATABLE :: dl_tmp1d … … 494 583 TYPE(TVAR) :: tl_lat1 495 584 496 TYPE(TFILE) :: tl_bathy1 497 TYPE(TFILE) :: tl_bathy0 498 499 TYPE(TMPP) :: tl_mppbathy1 500 TYPE(TMPP) :: tl_mppbathy0 585 TYPE(TMPP) :: tl_bathy1 586 TYPE(TMPP) :: tl_bathy0 501 587 502 588 TYPE(TDOM) :: tl_dom1 … … 510 596 IF( td_bdy%l_use )THEN 511 597 DO jl=1,td_bdy%i_nseg 512 513 !1- get boundary definition 598 ! get boundary definition 514 599 SELECT CASE(TRIM(td_bdy%c_card)) 515 600 CASE('north') … … 520 605 il_jmax1=td_bdy%t_seg(jl)%i_index 521 606 522 il_imin=1523 il_imax=il_imax1-il_imin1+1524 il_jmin=td_bdy%t_seg(jl)%i_width525 il_jmax=1526 527 607 CASE('south') 528 608 … … 532 612 il_jmax1=td_bdy%t_seg(jl)%i_index+(td_bdy%t_seg(jl)%i_width-1) 533 613 534 il_imin=1535 il_imax=il_imax1-il_imin1+1536 il_jmin=1537 il_jmax=td_bdy%t_seg(jl)%i_width538 539 614 CASE('east') 540 615 … … 544 619 il_jmax1=td_bdy%t_seg(jl)%i_last 545 620 546 il_imin=td_bdy%t_seg(jl)%i_width547 il_imax=1548 il_jmin=1549 il_jmax=il_jmax1-il_jmin1+1550 551 621 CASE('west') 552 622 … … 556 626 il_jmax1=td_bdy%t_seg(jl)%i_last 557 627 558 il_imin=1559 il_imax=td_bdy%t_seg(jl)%i_width560 il_jmin=1561 il_jmax=il_jmax1-il_jmin1+1562 563 628 END SELECT 564 629 565 !2 -read fine grid domain 566 tl_bathy1=td_bathy1 567 CALL iom_open(tl_bathy1) 568 569 !2-1 compute domain 630 ! -read fine grid domain 631 tl_bathy1=mpp_copy(td_bathy1) 632 633 ! compute domain 570 634 tl_dom1=dom_init( tl_bathy1, & 571 635 & il_imin1, il_imax1,& 572 & il_jmin1, il_jmax1 ) 573 574 !2-2 close file 575 CALL iom_close(tl_bathy1) 576 577 !2-3 read variables on domain (ugly way to do it, have to work on it) 578 !2-3-1 init mpp structure 579 tl_mppbathy1=mpp_init(tl_bathy1) 580 581 CALL file_clean(tl_bathy1) 582 583 !2-3-2 get processor to be used 584 CALL mpp_get_use( tl_mppbathy1, tl_dom1 ) 585 586 !2-3-3 open mpp files 587 CALL iom_mpp_open(tl_mppbathy1) 588 589 !2-3-4 read variable value on domain 590 tl_lon1=iom_mpp_read_var(tl_mppbathy1,'longitude',td_dom=tl_dom1) 591 tl_lat1=iom_mpp_read_var(tl_mppbathy1,'latitude' ,td_dom=tl_dom1) 592 593 !2-3-5 close mpp files 594 CALL iom_mpp_close(tl_mppbathy1) 595 596 !2-3-6 clean structure 597 CALL mpp_clean(tl_mppbathy1) 598 599 !3- get coarse grid indices 600 il_ind(:,:,:)=grid_get_coarse_index(td_bathy0, tl_lon1, tl_lat1, & 601 & id_rho=id_rho(:)) 602 603 il_imin0=il_ind(1,1,1) 604 il_imax0=il_ind(1,2,1) 605 606 il_jmin0=il_ind(2,1,1) 607 il_jmax0=il_ind(2,2,1) 608 609 il_offset(:,:)=il_ind(:,:,2) 610 611 !4- read coarse grid bathymetry on domain 612 tl_bathy0=td_bathy0 613 CALL iom_open(tl_bathy0) 614 615 !4-1 compute domain 636 & il_jmin1, il_jmax1,& 637 & TRIM(td_bdy%c_card)) 638 639 ! add extra band to fine grid domain (if possible) 640 ! to avoid dimension of one and so be able to compute offset 641 CALL dom_add_extra(tl_dom1, id_rho(jp_I), id_rho(jp_J)) 642 643 ! open mpp files over domain 644 CALL iom_dom_open(tl_bathy1, tl_dom1) 645 646 ! read variable value on domain 647 tl_lon1=iom_dom_read_var(tl_bathy1,'longitude',tl_dom1) 648 tl_lat1=iom_dom_read_var(tl_bathy1,'latitude' ,tl_dom1) 649 650 ! close mpp files 651 CALL iom_dom_close(tl_bathy1) 652 653 ! clean structure 654 CALL mpp_clean(tl_bathy1) 655 656 ! get coarse grid indices 657 il_ind(:,:)=grid_get_coarse_index(td_bathy0, tl_lon1, tl_lat1, & 658 & id_rho=id_rho(:)) 659 660 il_imin0=il_ind(1,1) 661 il_imax0=il_ind(1,2) 662 663 il_jmin0=il_ind(2,1) 664 il_jmax0=il_ind(2,2) 665 666 ! read coarse grid bathymetry on domain 667 tl_bathy0=mpp_copy(td_bathy0) 668 669 ! compute domain 616 670 tl_dom0=dom_init( tl_bathy0, & 617 671 & il_imin0, il_imax0,& 618 672 & il_jmin0, il_jmax0 ) 619 673 620 !4-2 close file 621 CALL iom_close(tl_bathy0) 622 623 !4-3 add extra band (if possible) to compute interpolation 674 il_offset(:,:)= grid_get_fine_offset(tl_bathy0, & 675 & il_imin0, il_jmin0,& 676 & il_imax0, il_jmax0,& 677 & tl_lon1%d_value(:,:,1,1), & 678 & tl_lat1%d_value(:,:,1,1), & 679 & id_rho=id_rho(:)) 680 681 ! clean 682 CALL var_clean(tl_lon1) 683 CALL var_clean(tl_lat1) 684 685 ! add extra band (if possible) to compute interpolation 624 686 CALL dom_add_extra(tl_dom0) 625 687 626 !4-4 read variables on domain (ugly way to do it, have to work on it) 627 !4-4-1 init mpp structure 628 tl_mppbathy0=mpp_init(tl_bathy0) 629 630 CALL file_clean(tl_bathy0) 631 632 !4-4-2 get processor to be used 633 CALL mpp_get_use( tl_mppbathy0, tl_dom0 ) 634 635 !4-4-3 open mpp files 636 CALL iom_mpp_open(tl_mppbathy0) 637 638 !4-4-4 read variable value on domain 639 tl_var0=iom_mpp_read_var(tl_mppbathy0,'Bathymetry',td_dom=tl_dom0) 640 641 !4-4-5 close mpp files 642 CALL iom_mpp_close(tl_mppbathy0) 643 644 !4-4-6 clean structure 645 CALL mpp_clean(tl_mppbathy0) 646 647 !5- interpolate variable 688 ! open mpp files over domain 689 CALL iom_dom_open(tl_bathy0, tl_dom0) 690 691 ! read variable value on domain 692 tl_var0=iom_dom_read_var(tl_bathy0,'Bathymetry',tl_dom0) 693 694 ! close mpp files 695 CALL iom_dom_close(tl_bathy0) 696 697 ! clean structure 698 CALL mpp_clean(tl_bathy0) 699 700 ! interpolate variable 648 701 CALL merge_bathy_interp( tl_var0, & 649 702 & id_rho(:), & 650 703 & il_offset(:,:) ) 651 704 652 ! 6-remove extraband added to domain705 ! remove extraband added to domain 653 706 CALL dom_del_extra( tl_var0, tl_dom0, id_rho(:) ) 654 707 655 ! 6-1remove extraband added to domain708 ! remove extraband added to domain 656 709 CALL dom_clean_extra( tl_dom0 ) 657 710 658 !7- fill refined table 659 !7-1 keep only useful point 660 ! interpolation could create more point than necessary 661 CALL boundary_clean_interp(tl_var0, td_bdy ) 662 663 ! use add request ???? 664 665 !7-2 fill refined table 711 ! remove extraband added to fine grid domain 712 CALL dom_del_extra( tl_var0, tl_dom1 ) 713 714 ! remove extraband added to fine grid domain 715 CALL dom_clean_extra( tl_dom1 ) 716 717 ! fill refined array 666 718 dd_refined( il_imin1:il_imax1, & 667 719 & il_jmin1:il_jmax1, & 668 720 & :,: )= tl_var0%d_value(:,:,:,:) 669 721 670 !8- compute weight function 722 ! clean 723 CALL var_clean(tl_var0) 724 725 ! compute weight function 671 726 ALLOCATE( dl_tmp1d(td_bdy%t_seg(jl)%i_width) ) 672 727 … … 678 733 679 734 ! compute weight on segment 680 dl_tmp1d(:)= 0.5 + 0.5*COS( (d g_pi*dl_tmp1d(:)) / &735 dl_tmp1d(:)= 0.5 + 0.5*COS( (dp_pi*dl_tmp1d(:)) / & 681 736 & (td_bdy%t_seg(jl)%i_width) ) 682 737 … … 694 749 695 750 ! compute weight on segment 696 dl_tmp1d(:)= 0.5 + 0.5*COS( (d g_pi*dl_tmp1d(:)) / &751 dl_tmp1d(:)= 0.5 + 0.5*COS( (dp_pi*dl_tmp1d(:)) / & 697 752 & (td_bdy%t_seg(jl)%i_width) ) 698 753 … … 710 765 711 766 ! compute weight on segment 712 dl_tmp1d(:)= 0.5 + 0.5*COS( (d g_pi*dl_tmp1d(:)) / &767 dl_tmp1d(:)= 0.5 + 0.5*COS( (dp_pi*dl_tmp1d(:)) / & 713 768 & (td_bdy%t_seg(jl)%i_width) ) 714 769 … … 726 781 727 782 ! compute weight on segment 728 dl_tmp1d(:)= 0.5 + 0.5*COS( (d g_pi*dl_tmp1d(:)) / &783 dl_tmp1d(:)= 0.5 + 0.5*COS( (dp_pi*dl_tmp1d(:)) / & 729 784 & (td_bdy%t_seg(jl)%i_width) ) 730 785 … … 740 795 DEALLOCATE( dl_tmp1d ) 741 796 742 ! 8-1 fill weight table797 ! fill weight array 743 798 ALLOCATE( dl_tmp2d( tl_dom1%t_dim(1)%i_len, & 744 799 & tl_dom1%t_dim(2)%i_len) ) … … 764 819 ENDIF 765 820 END SUBROUTINE merge_bathy_get_boundary 766 !> @endcode767 821 !------------------------------------------------------------------- 768 822 !> @brief 769 !> This subroutine 823 !> This subroutine interpolate variable. 770 824 !> 771 !> @details 825 !> @author J.Paul 826 !> @date November, 2013 - Initial Version 772 827 !> 773 !> @ author J.Paul774 !> - Nov, 2013- Initial Version775 !> 776 !> @param[in] 777 !> @ todo828 !> @param[inout] td_var variable structure 829 !> @param[in] id_rho array of refinment factor 830 !> @param[in] id_offset array of offset between fine and coarse grid 831 !> @param[in] id_iext i-direction size of extra bands (default=im_minext) 832 !> @param[in] id_jext j-direction size of extra bands (default=im_minext) 778 833 !------------------------------------------------------------------- 779 !> @code780 834 SUBROUTINE merge_bathy_interp( td_var, & 781 835 & id_rho, & … … 793 847 794 848 ! local variable 795 TYPE(TVAR) :: tl_var796 849 TYPE(TVAR) :: tl_mask 797 850 … … 803 856 ! loop indices 804 857 !---------------------------------------------------------------- 805 806 ! copy variable807 tl_var=td_var808 858 809 859 !WARNING: two extrabands are required for cubic interpolation … … 815 865 816 866 IF( il_iext < 2 .AND. td_var%c_interp(1) == 'cubic' )THEN 817 CALL logger_warn(" CREATE BATHY INTERP: at least extrapolation "//&867 CALL logger_warn("MERGE BATHY INTERP: at least extrapolation "//& 818 868 & "on two points are required with cubic interpolation ") 819 869 il_iext=2 … … 821 871 822 872 IF( il_jext < 2 .AND. td_var%c_interp(1) == 'cubic' )THEN 823 CALL logger_warn(" CREATE BATHY INTERP: at least extrapolation "//&873 CALL logger_warn("MERGE BATHY INTERP: at least extrapolation "//& 824 874 & "on two points are required with cubic interpolation ") 825 875 il_jext=2 826 876 ENDIF 827 877 828 ! 1-work on mask829 ! 1-1create mask830 ALLOCATE(bl_mask(t l_var%t_dim(1)%i_len, &831 & t l_var%t_dim(2)%i_len, &832 & t l_var%t_dim(3)%i_len, &833 & t l_var%t_dim(4)%i_len) )878 ! work on mask 879 ! create mask 880 ALLOCATE(bl_mask(td_var%t_dim(1)%i_len, & 881 & td_var%t_dim(2)%i_len, & 882 & td_var%t_dim(3)%i_len, & 883 & td_var%t_dim(4)%i_len) ) 834 884 835 885 bl_mask(:,:,:,:)=1 836 WHERE(t l_var%d_value(:,:,:,:)==tl_var%d_fill) bl_mask(:,:,:,:)=0837 838 SELECT CASE(TRIM(t l_var%c_point))886 WHERE(td_var%d_value(:,:,:,:)==td_var%d_fill) bl_mask(:,:,:,:)=0 887 888 SELECT CASE(TRIM(td_var%c_point)) 839 889 CASE DEFAULT ! 'T' 840 tl_mask=var_init('tmask', bl_mask(:,:,:,:), td_dim=tl_var%t_dim(:))841 CASE('U')842 tl_mask=var_init('umask', bl_mask(:,:,:,:), td_dim=tl_var%t_dim(:))843 CASE('V')844 tl_mask=var_init('vmask', bl_mask(:,:,:,:), td_dim=tl_var%t_dim(:))845 CASE('F')846 tl_mask=var_init('fmask', bl_mask(:,:,:,:), td_dim=tl_var%t_dim(:))890 tl_mask=var_init('tmask',bl_mask(:,:,:,:),td_dim=td_var%t_dim(:),& 891 & id_ew=td_var%i_ew ) 892 CASE('U','V','F') 893 CALL logger_fatal("MERGE BATHY INTERP: can not computed "//& 894 & "interpolation on "//TRIM(td_var%c_point)//& 895 & " grid point (variable "//TRIM(td_var%c_name)//& 896 & "). check namelist.") 847 897 END SELECT 848 898 849 899 DEALLOCATE(bl_mask) 850 900 851 ! 1-2interpolate mask901 ! interpolate mask 852 902 CALL interp_fill_value( tl_mask, id_rho(:), & 853 903 & id_offset=id_offset(:,:) ) 854 904 855 ! 2-work on variable856 ! 2-0add extraband857 CALL extrap_add_extrabands(t l_var, il_iext, il_iext)858 859 ! 2-1extrapolate variable860 CALL extrap_fill_value( t l_var, id_offset=id_offset(:,:), &905 ! work on variable 906 ! add extraband 907 CALL extrap_add_extrabands(td_var, il_iext, il_iext) 908 909 ! extrapolate variable 910 CALL extrap_fill_value( td_var, id_offset=id_offset(:,:), & 861 911 & id_rho=id_rho(:), & 862 912 & id_iext=il_iext, id_jext=il_jext ) 863 913 864 ! 2-2interpolate Bathymetry865 CALL interp_fill_value( t l_var, id_rho(:), &914 ! interpolate Bathymetry 915 CALL interp_fill_value( td_var, id_rho(:), & 866 916 & id_offset=id_offset(:,:) ) 867 917 868 ! 2-3remove extraband869 CALL extrap_del_extrabands(t l_var, il_iext*id_rho(jp_I), il_jext*id_rho(jp_J))870 871 ! 2-4keep original mask918 ! remove extraband 919 CALL extrap_del_extrabands(td_var, il_iext*id_rho(jp_I), il_jext*id_rho(jp_J)) 920 921 ! keep original mask 872 922 WHERE( tl_mask%d_value(:,:,:,:) == 0 ) 873 t l_var%d_value(:,:,:,:)=tl_var%d_fill923 td_var%d_value(:,:,:,:)=td_var%d_fill 874 924 END WHERE 875 925 876 !3- save result877 td_var=tl_var878 879 ! clean variable structure880 CALL var_clean(tl_var)881 882 926 END SUBROUTINE merge_bathy_interp 883 !> @endcode884 927 END PROGRAM merge_bathy
Note: See TracChangeset
for help on using the changeset viewer.
