- Timestamp:
- 2015-07-16T13:55:15+02:00 (9 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/TOOLS/SIREN/src/create_coord.f90
r4213 r5602 7 7 ! 8 8 ! DESCRIPTION: 9 !> @file 9 10 !> @brief 10 !> This program create coordinate file.11 !> This program create fine grid coordinate file. 11 12 !> 12 13 !> @details 13 !> Variables are extracted from the input coordinates coarse grid, 14 !> and interpolated to create fine coordinates files. 15 !> 16 !> @author 17 !> J.Paul 14 !> @section sec1 method 15 !> All variables from the input coordinates coarse grid file, are extracted 16 !> and interpolated to create fine grid coordinates files.<br/> 17 !> @note 18 !> interpolation method could be different for each variable. 19 !> 20 !> @section sec2 how to 21 !> to create fine grid coordinates files:<br/> 22 !> @code{.sh} 23 !> ./SIREN/bin/create_coord create_coord.nam 24 !> @endcode 25 !> 26 !> create_coord.nam comprise 6 namelists:<br/> 27 !> - logger namelist (namlog) 28 !> - config namelist (namcfg) 29 !> - coarse grid namelist (namcrs) 30 !> - variable namelist (namvar) 31 !> - nesting namelist (namnst) 32 !> - output namelist (namout) 33 !> 34 !> @note 35 !> All namelists have to be in file create_coord.nam, 36 !> however variables of those namelists are all optional. 37 !> 38 !> * _logger namelist (namlog)_:<br/> 39 !> - cn_logfile : log filename 40 !> - cn_verbosity : verbosity ('trace','debug','info', 41 !> 'warning','error','fatal') 42 !> - in_maxerror : maximum number of error allowed 43 !> 44 !> * _config namelist (namcfg)_:<br/> 45 !> - cn_varcfg : variable configuration file 46 !> (see ./SIREN/cfg/variable.cfg) 47 !> 48 !> * _coarse grid namelist (namcrs)_:<br/> 49 !> - cn_coord0 : coordinate file 50 !> - in_perio0 : NEMO periodicity index (see Model Boundary Condition in 51 !> [NEMO documentation](http://www.nemo-ocean.eu/About-NEMO/Reference-manuals)) 52 !> 53 !> * _variable namelist (namvar)_:<br/> 54 !> - cn_varinfo : list of variable and extra information about request(s) 55 !> to be used.<br/> 56 !> each elements of *cn_varinfo* is a string character.<br/> 57 !> it is composed of the variable name follow by ':', 58 !> then request(s) to be used on this variable.<br/> 59 !> request could be: 60 !> - interpolation method 61 !> - extrapolation method 62 !> - filter method 63 !> 64 !> requests must be separated by ';' .<br/> 65 !> order of requests does not matter.<br/> 66 !> 67 !> informations about available method could be find in @ref interp, 68 !> @ref extrap and @ref filter modules.<br/> 69 !> 70 !> Example: 'votemper: linear; hann(2,3); dist_weight', 71 !> 'vosaline: cubic'<br/> 72 !> @note 73 !> If you do not specify a method which is required, 74 !> default one is applied. 75 !> 76 !> * _nesting namelist (namnst)_:<br/> 77 !> - in_imin0 : i-direction lower left point indice 78 !> of coarse grid subdomain to be used 79 !> - in_imax0 : i-direction upper right point indice 80 !> of coarse grid subdomain to be used 81 !> - in_jmin0 : j-direction lower left point indice 82 !> of coarse grid subdomain to be used 83 !> - in_jmax0 : j-direction upper right point indice 84 !> of coarse grid subdomain to be used 85 !> - in_rhoi : refinement factor in i-direction 86 !> - in_rhoj : refinement factor in j-direction<br/> 87 !> 88 !> \image html grid_zoom_40.png 89 !> \image latex grid_zoom_40.png 90 !> 91 !> * _output namelist (namout)_: 92 !> - cn_fileout : output coordinate file 93 !> 94 !> @author J.Paul 18 95 ! REVISION HISTORY: 19 !> @date Nov, 2013 - Initial Version 20 ! 96 !> @date November, 2013 - Initial Version 97 !> @date September, 2014 98 !> - add header for user 99 !> - compute offset considering grid point 100 !> - add global attributes in output file 101 !> 21 102 !> @note Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 22 !>23 !> @todo24 !> - add extrapolation (case coordin with mask)25 !> - add extraction from a grid at fine resolution26 103 !---------------------------------------------------------------------- 27 !> @code28 104 PROGRAM create_coord 29 105 30 ! USE netcdf ! nf90 library31 106 USE global ! global variable 32 107 USE kind ! F90 kind parameter … … 39 114 USE file ! file manager 40 115 USE iom ! I/O manager 41 USE dom ! domain manager42 116 USE grid ! grid manager 43 117 USE extrap ! extrapolation manager … … 45 119 USE filter ! filter manager 46 120 USE mpp ! MPP manager 121 USE dom ! domain manager 47 122 USE iom_mpp ! MPP I/O manager 123 USE iom_dom ! DOM I/O manager 48 124 49 125 IMPLICIT NONE … … 56 132 INTEGER(i4) :: il_status 57 133 INTEGER(i4) :: il_fileid 134 INTEGER(i4) :: il_attid 135 INTEGER(i4) :: il_ind 58 136 INTEGER(i4) :: il_nvar 59 ! INTEGER(i4) , DIMENSION(:,:,:,:) , ALLOCATABLE :: il_value 137 INTEGER(i4) :: il_ew 60 138 INTEGER(i4) , DIMENSION(ip_maxdim) :: il_rho 61 139 INTEGER(i4) , DIMENSION(2,2,ip_npoint) :: il_offset 62 140 63 141 LOGICAL :: ll_exist … … 71 149 TYPE(TDIM) , DIMENSION(ip_maxdim) :: tl_dim 72 150 73 TYPE(T FILE):: tl_coord0151 TYPE(TMPP) :: tl_coord0 74 152 TYPE(TFILE) :: tl_fileout 75 153 76 TYPE(TMPP) :: tl_mppcoordin 154 ! check 155 ! INTEGER(i4) :: il_imin0 156 ! INTEGER(i4) :: il_imax0 157 ! INTEGER(i4) :: il_jmin0 158 ! INTEGER(i4) :: il_jmax0 159 ! INTEGER(i4) , DIMENSION(2,2) :: il_ind2 160 ! TYPE(TMPP) :: tl_mppout 77 161 78 162 ! loop indices 79 163 INTEGER(i4) :: ji 164 INTEGER(i4) :: jj 80 165 81 166 ! namelist variable 82 167 CHARACTER(LEN=lc) :: cn_logfile = 'create_coord.log' 83 168 CHARACTER(LEN=lc) :: cn_verbosity = 'warning' 169 INTEGER(i4) :: in_maxerror = 5 84 170 85 171 CHARACTER(LEN=lc) :: cn_coord0 = '' 86 172 INTEGER(i4) :: in_perio0 = -1 87 173 88 CHARACTER(LEN=lc) :: cn_varcfg = ' variable.cfg'89 90 CHARACTER(LEN=lc), DIMENSION(i g_maxvar) :: cn_varinfo = ''174 CHARACTER(LEN=lc) :: cn_varcfg = '../cfg/variable.cfg' 175 176 CHARACTER(LEN=lc), DIMENSION(ip_maxvar) :: cn_varinfo = '' 91 177 92 178 INTEGER(i4) :: in_imin0 = 0 … … 100 186 !------------------------------------------------------------------- 101 187 102 NAMELIST /namlog/ & !< logger namelist 103 & cn_logfile, & !< log file 104 & cn_verbosity !< logger verbosity 105 106 NAMELIST /namcfg/ & !< config namelist 188 NAMELIST /namlog/ & ! logger namelist 189 & cn_logfile, & !< logger file name 190 & cn_verbosity, & !< logger verbosity 191 & in_maxerror !< logger maximum error 192 193 NAMELIST /namcfg/ & ! config namelist 107 194 & cn_varcfg !< variable configuration file 108 195 109 NAMELIST /namcrs/ & ! coarse grid namelist196 NAMELIST /namcrs/ & ! coarse grid namelist 110 197 & cn_coord0 , & !< coordinate file 111 198 & in_perio0 !< periodicity index 112 199 113 NAMELIST /namvar/ & ! namvar200 NAMELIST /namvar/ & ! variable namelist 114 201 & cn_varinfo !< list of variable and extra information about 115 202 !< interpolation, extrapolation or filter method to be used. 116 !< (ex: 'votemper /linear/hann/dist_weight','vosaline/cubic' )203 !< (ex: 'votemper:linear,hann,dist_weight','vosaline:cubic' ) 117 204 118 NAMELIST /namnst/ & ! <nesting namelist205 NAMELIST /namnst/ & ! nesting namelist 119 206 & in_imin0, & !< i-direction lower left point indice 120 207 & in_imax0, & !< i-direction upper right point indice … … 124 211 & in_rhoj !< refinement factor in j-direction 125 212 126 NAMELIST /namout/ & ! <output namelist127 & cn_fileout !< fine grid coordinate file213 NAMELIST /namout/ & ! output namelist 214 & cn_fileout !< fine grid coordinate file 128 215 !------------------------------------------------------------------- 129 216 130 ! 1-namelist131 ! 1-1get namelist217 ! namelist 218 ! get namelist 132 219 il_narg=COMMAND_ARGUMENT_COUNT() !f03 intrinsec 133 220 IF( il_narg/=1 )THEN … … 138 225 ENDIF 139 226 140 ! 1-2read namelist227 ! read namelist 141 228 INQUIRE(FILE=TRIM(cl_namelist), EXIST=ll_exist) 142 229 IF( ll_exist )THEN 143 230 144 231 il_fileid=fct_getunit() 145 232 … … 157 244 158 245 READ( il_fileid, NML = namlog ) 159 ! 1-2-1 define logfile160 CALL logger_open(TRIM(cn_logfile),TRIM(cn_verbosity) )246 ! define logger file 247 CALL logger_open(TRIM(cn_logfile),TRIM(cn_verbosity),in_maxerror) 161 248 CALL logger_header() 162 249 163 250 READ( il_fileid, NML = namcfg ) 164 ! 1-2-2get variable extra information on configuration file251 ! get variable extra information on configuration file 165 252 CALL var_def_extra(TRIM(cn_varcfg)) 166 253 167 254 READ( il_fileid, NML = namcrs ) 168 255 READ( il_fileid, NML = namvar ) 169 ! 1-2-3add user change in extra information256 ! add user change in extra information 170 257 CALL var_chg_extra( cn_varinfo ) 171 258 … … 182 269 183 270 PRINT *,"ERROR in create_coord: can't find "//TRIM(cl_namelist) 184 185 ENDIF 186 187 !2- open files 271 STOP 272 273 ENDIF 274 275 ! open files 188 276 IF( cn_coord0 /= '' )THEN 189 tl_coord0= file_init(TRIM(cn_coord0),id_perio=in_perio0)190 CALL iom_open(tl_coord0)277 tl_coord0=mpp_init( file_init(TRIM(cn_coord0)), id_perio=in_perio0) 278 CALL grid_get_info(tl_coord0) 191 279 ELSE 192 280 CALL logger_fatal("CREATE COORD: no coarse grid coordinate found. "//& 193 281 & "check namelist") 194 282 ENDIF 195 283 196 ! 3-check197 ! 3-1check output file do not already exist284 ! check 285 ! check output file do not already exist 198 286 INQUIRE(FILE=TRIM(cn_fileout), EXIST=ll_exist) 199 287 IF( ll_exist )THEN … … 202 290 ENDIF 203 291 204 ! 3-2 check namelist205 IF( in_imin0 < 1 .OR. in_imax0 < 1 .OR. in_jmin0 < 1 .OR. in_jmax0 < 1)THEN206 CALL logger_ error("CREATE COORD: invalid point indice."//&292 ! check nesting parameters 293 IF( in_imin0 < 0 .OR. in_imax0 < 0 .OR. in_jmin0 < 0 .OR. in_jmax0 < 0)THEN 294 CALL logger_fatal("CREATE COORD: invalid points indices."//& 207 295 & " check namelist "//TRIM(cl_namelist)) 208 296 ENDIF … … 215 303 il_rho(jp_I)=in_rhoi 216 304 il_rho(jp_J)=in_rhoj 217 ENDIF 218 219 !3-3 check domain validity 305 306 il_offset(:,:,:)=create_coord_get_offset(il_rho(:)) 307 308 ENDIF 309 310 ! check domain validity 220 311 CALL grid_check_dom(tl_coord0, in_imin0, in_imax0, in_jmin0, in_jmax0 ) 221 312 222 ! 4-compute domain313 ! compute domain 223 314 tl_dom=dom_init( tl_coord0, & 224 315 & in_imin0, in_imax0,& 225 316 & in_jmin0, in_jmax0 ) 226 317 227 ! close file 228 CALL iom_close(tl_coord0) 229 230 !4-1 add extra band (if possible) to compute interpolation 318 ! add extra band (if need be) to compute interpolation 231 319 CALL dom_add_extra(tl_dom) 232 320 233 !5- read variables on domain (ugly way to do it, have to work on it) 234 !5-1 init mpp structure 235 tl_mppcoordin=mpp_init(tl_coord0) 236 237 CALL file_clean(tl_coord0) 238 239 !5-2 get processor to be used 240 CALL mpp_get_use( tl_mppcoordin, tl_dom ) 241 242 !5-3 open mpp files 243 CALL iom_mpp_open(tl_mppcoordin) 244 245 !5-4 fill variable value on domain 246 CALL iom_mpp_fill_var(tl_mppcoordin, tl_dom) 247 248 !5-5 close mpp files 249 CALL iom_mpp_close(tl_mppcoordin) 250 251 il_nvar=tl_mppcoordin%t_proc(1)%i_nvar 321 ! open mpp files 322 CALL iom_dom_open(tl_coord0, tl_dom) 323 324 il_nvar=tl_coord0%t_proc(1)%i_nvar 252 325 ALLOCATE( tl_var(il_nvar) ) 253 326 DO ji=1,il_nvar 254 327 255 tl_var(ji)=tl_mppcoordin%t_proc(1)%t_var(ji) 256 !7- interpolate variables 257 CALL create_coord_interp( tl_var(ji), il_rho(:) ) 258 259 !6- remove extraband added to domain 260 CALL dom_del_extra( tl_var(ji), tl_dom, il_rho(:) ) 261 262 !7- add ghost cell 263 CALL grid_add_ghost(tl_var(ji),tl_dom%i_ighost,tl_dom%i_jghost) 264 265 !8- filter 328 tl_var(ji)=iom_dom_read_var(tl_coord0, & 329 & TRIM(tl_coord0%t_proc(1)%t_var(ji)%c_name),& 330 & tl_dom) 331 332 SELECT CASE(TRIM(tl_var(ji)%c_point)) 333 CASE('T') 334 jj=jp_T 335 CASE('U') 336 jj=jp_U 337 CASE('V') 338 jj=jp_V 339 CASE('F') 340 jj=jp_F 341 END SELECT 342 343 ! interpolate variables 344 CALL create_coord_interp( tl_var(ji), il_rho(:), & 345 & il_offset(:,:,jj) ) 346 347 ! remove extraband added to domain 348 CALL dom_del_extra( tl_var(ji), tl_dom, il_rho(:), .true. ) 349 350 ! do not add ghost cell. 351 ! ghost cell already replace by value for coordinates 352 ! CALL grid_add_ghost(tl_var(ji),tl_dom%i_ghost(:,:)) 353 354 ! filter 266 355 CALL filter_fill_value(tl_var(ji)) 267 356 268 357 ENDDO 269 358 270 ! 9- clean271 DO ji=1,il_nvar272 CALL var_clean(tl_mppcoordin%t_proc(1)%t_var(ji)) 273 ENDDO274 CALL mpp_clean(tl_ mppcoordin)275 276 ! 10-create file359 ! close mpp files 360 CALL iom_dom_close(tl_coord0) 361 362 ! clean 363 CALL mpp_clean(tl_coord0) 364 365 ! create file 277 366 tl_fileout=file_init(TRIM(cn_fileout)) 278 367 279 ! 10-1add dimension368 ! add dimension 280 369 ! save biggest dimension 281 370 tl_dim(:)=var_max_dim(tl_var(:)) … … 285 374 ENDDO 286 375 287 !10-2 add variables 288 376 ! add variables 289 377 DO ji=1,il_nvar 290 378 CALL file_add_var(tl_fileout, tl_var(ji)) 291 379 ENDDO 292 380 293 !10-3 add some attribute 381 ! recompute some attribute 382 383 ! add some attribute 294 384 tl_att=att_init("Created_by","SIREN create_coord") 295 385 CALL file_add_att(tl_fileout, tl_att) … … 299 389 CALL file_add_att(tl_fileout, tl_att) 300 390 301 tl_att=att_init("s ource_file",TRIM(fct_basename(cn_coord0)))391 tl_att=att_init("src_file",TRIM(fct_basename(cn_coord0))) 302 392 CALL file_add_att(tl_fileout, tl_att) 303 393 304 tl_att=att_init("s ource_i-indices",(/in_imin0,in_imax0/))394 tl_att=att_init("src_i_indices",(/in_imin0,in_imax0/)) 305 395 CALL file_add_att(tl_fileout, tl_att) 306 tl_att=att_init("source_j-indices",(/in_jmin0,in_jmax0/)) 307 CALL file_add_att(tl_fileout, tl_att) 308 309 !10-4 create file 396 tl_att=att_init("src_j_indices",(/in_jmin0,in_jmax0/)) 397 CALL file_add_att(tl_fileout, tl_att) 398 IF( .NOT. ALL(il_rho(:)==1) )THEN 399 tl_att=att_init("refinment_factor",(/il_rho(jp_I),il_rho(jp_J)/)) 400 CALL file_add_att(tl_fileout, tl_att) 401 ENDIF 402 403 ! add attribute periodicity 404 il_attid=0 405 IF( ASSOCIATED(tl_fileout%t_att) )THEN 406 il_attid=att_get_index(tl_fileout%t_att(:),'periodicity') 407 ENDIF 408 IF( tl_dom%i_perio >= 0 .AND. il_attid == 0 )THEN 409 tl_att=att_init('periodicity',tl_dom%i_perio) 410 CALL file_add_att(tl_fileout,tl_att) 411 ENDIF 412 413 ! add attribute east west overlap 414 il_attid=0 415 IF( ASSOCIATED(tl_fileout%t_att) )THEN 416 il_attid=att_get_index(tl_fileout%t_att(:),'ew_overlap') 417 ENDIF 418 IF( il_attid == 0 )THEN 419 il_ind=var_get_index(tl_fileout%t_var(:),'longitude') 420 il_ew=grid_get_ew_overlap(tl_fileout%t_var(il_ind)) 421 IF( il_ew >= 0 )THEN 422 tl_att=att_init('ew_overlap',il_ew) 423 CALL file_add_att(tl_fileout,tl_att) 424 ENDIF 425 ENDIF 426 427 ! create file 310 428 CALL iom_create(tl_fileout) 311 429 312 ! 10-5write file430 ! write file 313 431 CALL iom_write_file(tl_fileout) 314 432 315 ! 10-6close file433 ! close file 316 434 CALL iom_close(tl_fileout) 317 435 318 !11- clean 319 DO ji=1,il_nvar 320 CALL var_clean(tl_var(ji)) 321 ENDDO 436 ! clean 437 CALL att_clean(tl_att) 438 CALL var_clean(tl_var(:)) 439 DEALLOCATE( tl_var) 440 322 441 CALL file_clean(tl_fileout) 323 442 324 DEALLOCATE( tl_var) 443 ! ! check domain 444 ! tl_coord0=mpp_init( file_init(TRIM(cn_coord0)), id_perio=in_perio0) 445 ! tl_mppout=mpp_init( file_init(TRIM(cn_fileout)) ) 446 ! CALL grid_get_info(tl_coord0) 447 ! CALL iom_mpp_open(tl_mppout) 448 ! 449 ! il_ind2(:,:)=grid_get_coarse_index( tl_coord0, tl_mppout, & 450 ! & id_rho=il_rho(:) ) 451 ! 452 ! il_imin0=il_ind2(1,1) ; il_imax0=il_ind2(1,2) 453 ! il_jmin0=il_ind2(2,1) ; il_jmax0=il_ind2(2,2) 454 ! 455 ! IF( il_imin0 /= in_imin0 .OR. & 456 ! & il_imax0 /= in_imax0 .OR. & 457 ! & il_jmin0 /= in_jmin0 .OR. & 458 ! & il_jmax0 /= in_jmax0 )THEN 459 ! CALL logger_debug("CREATE COORD: output indices ("//& 460 ! & TRIM(fct_str(il_imin0))//","//& 461 ! & TRIM(fct_str(il_imax0))//") ("//& 462 ! & TRIM(fct_str(il_jmin0))//","//& 463 ! & TRIM(fct_str(il_jmax0))//")" ) 464 ! CALL logger_debug("CREATE COORD: input indices ("//& 465 ! & TRIM(fct_str(in_imin0))//","//& 466 ! & TRIM(fct_str(in_imax0))//") ("//& 467 ! & TRIM(fct_str(in_jmin0))//","//& 468 ! & TRIM(fct_str(in_jmax0))//")" ) 469 ! CALL logger_fatal("CREATE COORD: output domain not confrom "//& 470 ! & "with input indices") 471 ! ENDIF 472 ! 473 ! CALL iom_mpp_close(tl_coord0) 474 ! CALL iom_mpp_close(tl_mppout) 325 475 326 476 ! close log file 327 477 CALL logger_footer() 328 CALL logger_close() 329 330 !> @endcode 478 CALL logger_close() 479 331 480 CONTAINS 332 481 !------------------------------------------------------------------- 333 482 !> @brief 334 !> This subroutine 483 !> This function compute offset over Arakawa grid points, 484 !> given refinement factor. 485 !> 486 !> @author J.Paul 487 !> @date August, 2014 - Initial Version 488 !> 489 !> @param[in] id_rho array of refinement factor 490 !> @return array of offset 491 !------------------------------------------------------------------- 492 FUNCTION create_coord_get_offset( id_rho ) 493 IMPLICIT NONE 494 ! Argument 495 INTEGER(i4), DIMENSION(:), INTENT(IN) :: id_rho 496 497 ! function 498 INTEGER(i4), DIMENSION(2,2,ip_npoint) :: create_coord_get_offset 499 ! local variable 500 ! loop indices 501 !---------------------------------------------------------------- 502 503 ! case 'T' 504 create_coord_get_offset(jp_I,:,jp_T)=FLOOR(REAL(id_rho(jp_I)-1,dp)*0.5) 505 create_coord_get_offset(jp_J,:,jp_T)=FLOOR(REAL(id_rho(jp_J)-1,dp)*0.5) 506 ! case 'U' 507 create_coord_get_offset(jp_I,1,jp_U)=0 508 create_coord_get_offset(jp_I,2,jp_U)=id_rho(jp_I)-1 509 create_coord_get_offset(jp_J,:,jp_U)=FLOOR(REAL(id_rho(jp_J)-1,dp)*0.5) 510 ! case 'V' 511 create_coord_get_offset(jp_I,:,jp_V)=FLOOR(REAL(id_rho(jp_I)-1,dp)*0.5) 512 create_coord_get_offset(jp_J,1,jp_V)=0 513 create_coord_get_offset(jp_J,2,jp_V)=id_rho(jp_J)-1 514 ! case 'F' 515 create_coord_get_offset(jp_I,1,jp_F)=0 516 create_coord_get_offset(jp_I,2,jp_F)=id_rho(jp_I)-1 517 create_coord_get_offset(jp_J,1,jp_F)=0 518 create_coord_get_offset(jp_J,2,jp_F)=id_rho(jp_J)-1 519 520 521 END FUNCTION create_coord_get_offset 522 !------------------------------------------------------------------- 523 !> @brief 524 !> This subroutine interpolate variable, given refinment factor. 335 525 !> 336 526 !> @details 527 !> Optionaly, you could specify number of points 528 !> to be extrapolated in i- and j-direction.<br/> 529 !> variable mask is first computed (using _FillValue) and interpolated.<br/> 530 !> variable is then extrapolated, and interpolated.<br/> 531 !> Finally interpolated mask is applied on refined variable. 337 532 !> 338 533 !> @author J.Paul 339 !> - Nov, 2013- Initial Version534 !> @date November, 2013 - Initial Version 340 535 !> 341 !> @param[in] 342 !> @todo 536 !> @param[inout] td_var variable strcuture 537 !> @param[in] id_rho array of refinement factor 538 !> @param[in] id_offset offset between fine grid and coarse grid 539 !> @param[in] id_iext number of points to be extrapolated in i-direction 540 !> @param[in] id_jext number of points to be extrapolated in j-direction 343 541 !------------------------------------------------------------------- 344 !> @code345 542 SUBROUTINE create_coord_interp( td_var, & 346 543 & id_rho, & 544 & id_offset, & 347 545 & id_iext, id_jext) 348 546 … … 350 548 351 549 ! Argument 352 TYPE(TVAR) , INTENT(INOUT) :: td_var 353 INTEGER(i4), DIMENSION(:), INTENT(IN ) :: id_rho 354 INTEGER(i4), INTENT(IN ), OPTIONAL :: id_iext 355 INTEGER(i4), INTENT(IN ), OPTIONAL :: id_jext 550 TYPE(TVAR) , INTENT(INOUT) :: td_var 551 INTEGER(i4), DIMENSION(:) , INTENT(IN ) :: id_rho 552 INTEGER(i4), DIMENSION(:,:), INTENT(IN) :: id_offset 553 INTEGER(i4), INTENT(IN ), OPTIONAL :: id_iext 554 INTEGER(i4), INTENT(IN ), OPTIONAL :: id_jext 356 555 357 556 ! local variable 358 557 TYPE(TVAR) :: tl_mask 359 TYPE(TVAR) :: tl_var360 558 361 559 INTEGER(i1), DIMENSION(:,:,:,:), ALLOCATABLE :: bl_mask 362 363 INTEGER(i4), DIMENSION(2,2) :: il_offset364 560 365 561 INTEGER(i4) :: il_iext … … 369 565 !---------------------------------------------------------------- 370 566 371 ! copy variable 372 tl_var=td_var 567 IF( ANY(SHAPE(id_offset(:,:)) /= 2) )THEN 568 CALL logger_error("CREATE COORD INTERP: invalid dimension of "//& 569 & "offset array") 570 ENDIF 373 571 374 572 !WARNING: two extrabands are required for cubic interpolation … … 391 589 ENDIF 392 590 393 !1- work on mask 394 !1-1 create mask 395 ALLOCATE(bl_mask(tl_var%t_dim(1)%i_len, & 396 & tl_var%t_dim(2)%i_len, & 397 & tl_var%t_dim(3)%i_len, & 398 & tl_var%t_dim(4)%i_len) ) 399 400 bl_mask(:,:,:,:)=1 401 WHERE(tl_var%d_value(:,:,:,:)==tl_var%d_fill) bl_mask(:,:,:,:)=0 402 403 SELECT CASE(TRIM(tl_var%c_point)) 404 CASE DEFAULT ! 'T' 405 tl_mask=var_init('tmask', bl_mask(:,:,:,:), td_dim=tl_var%t_dim(:)) 406 CASE('U') 407 tl_mask=var_init('umask', bl_mask(:,:,:,:), td_dim=tl_var%t_dim(:)) 408 CASE('V') 409 tl_mask=var_init('vmask', bl_mask(:,:,:,:), td_dim=tl_var%t_dim(:)) 410 CASE('F') 411 tl_mask=var_init('fmask', bl_mask(:,:,:,:), td_dim=tl_var%t_dim(:)) 412 END SELECT 413 414 DEALLOCATE(bl_mask) 415 416 !1-2 interpolate mask 417 il_offset(:,:)=1 418 CALL interp_fill_value( tl_mask, id_rho(:), & 419 & id_offset=il_offset(:,:) ) 420 421 !2- work on variable 422 !2-0 add extraband 423 CALL extrap_add_extrabands(tl_var, il_iext, il_jext) 424 425 !2-1 extrapolate variable 426 CALL extrap_fill_value( tl_var, id_iext=il_iext, id_jext=il_jext ) 427 428 !2-2 interpolate variable 429 il_offset(:,:)=1 430 CALL interp_fill_value( tl_var, id_rho(:), & 431 & id_offset=il_offset(:,:)) 432 433 !2-3 remove extraband 434 CALL extrap_del_extrabands(tl_var, il_iext*id_rho(jp_I), il_jext*id_rho(jp_J)) 435 436 !3- keep original mask 437 WHERE( tl_mask%d_value(:,:,:,:) == 0 ) 438 tl_var%d_value(:,:,:,:)=tl_var%d_fill 439 END WHERE 440 441 !4- save result 442 td_var=tl_var 591 IF( ANY(id_rho(:)>1) )THEN 592 ! work on mask 593 ! create mask 594 ALLOCATE(bl_mask(td_var%t_dim(1)%i_len, & 595 & td_var%t_dim(2)%i_len, & 596 & td_var%t_dim(3)%i_len, & 597 & td_var%t_dim(4)%i_len) ) 598 599 bl_mask(:,:,:,:)=1 600 WHERE(td_var%d_value(:,:,:,:)==td_var%d_fill) bl_mask(:,:,:,:)=0 601 602 SELECT CASE(TRIM(td_var%c_point)) 603 CASE DEFAULT ! 'T' 604 tl_mask=var_init('tmask',bl_mask(:,:,:,:),td_dim=td_var%t_dim(:),& 605 & id_ew=td_var%i_ew ) 606 CASE('U') 607 tl_mask=var_init('umask',bl_mask(:,:,:,:),td_dim=td_var%t_dim(:),& 608 & id_ew=td_var%i_ew ) 609 CASE('V') 610 tl_mask=var_init('vmask',bl_mask(:,:,:,:),td_dim=td_var%t_dim(:),& 611 & id_ew=td_var%i_ew ) 612 CASE('F') 613 tl_mask=var_init('fmask',bl_mask(:,:,:,:),td_dim=td_var%t_dim(:),& 614 & id_ew=td_var%i_ew ) 615 END SELECT 616 617 DEALLOCATE(bl_mask) 618 619 ! interpolate mask 620 CALL interp_fill_value( tl_mask, id_rho(:), & 621 & id_offset=id_offset(:,:) ) 622 623 ! work on variable 624 ! add extraband 625 CALL extrap_add_extrabands(td_var, il_iext, il_jext) 626 627 ! extrapolate variable 628 CALL extrap_fill_value( td_var, id_iext=il_iext, id_jext=il_jext ) 629 630 ! interpolate variable 631 CALL interp_fill_value( td_var, id_rho(:), & 632 & id_offset=id_offset(:,:)) 633 634 ! remove extraband 635 CALL extrap_del_extrabands(td_var, il_iext*id_rho(jp_I), il_jext*id_rho(jp_J)) 636 637 ! keep original mask 638 WHERE( tl_mask%d_value(:,:,:,:) == 0 ) 639 td_var%d_value(:,:,:,:)=td_var%d_fill 640 END WHERE 641 ENDIF 443 642 444 643 ! clean variable structure 445 644 CALL var_clean(tl_mask) 446 CALL var_clean(tl_var)447 645 448 646 END SUBROUTINE create_coord_interp 449 !> @endcode450 647 END PROGRAM create_coord
Note: See TracChangeset
for help on using the changeset viewer.
