Changeset 3841 for palm/trunk/UTIL
- Timestamp:
- Mar 29, 2019 11:38:44 AM (6 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
palm/trunk/UTIL/surface_output_processing/surface_output_to_vtk.f90
r3805 r3841 25 25 ! ----------------- 26 26 ! $Id$ 27 ! Add suffix to VTK output to indicate average data 28 ! 29 ! 3805 2019-03-20 15:26:35Z raasch 27 30 ! output format adjusted 28 31 ! … … 143 146 !> Output is distinguished between instantaneous and time-averaged data. 144 147 !------------------------------------------------------------------------------! 145 PROGRAM surface_output_to_vtk 148 PROGRAM surface_output_to_vtk_netcdf 146 149 147 150 USE, INTRINSIC :: ISO_C_BINDING 151 152 #if defined( __netcdf ) 153 USE NETCDF 154 #endif 148 155 149 156 USE posix_interface, & … … 171 178 INTEGER, PARAMETER :: OFFSET_KIND = C_SIZE_T !< unsigned integer for the C-interface 172 179 173 INTEGER(iwp) :: cycle_number 180 INTEGER(iwp) :: cycle_number = 0 !< cycle number 174 181 INTEGER(iwp) :: f !< running index over all binary files 175 182 INTEGER(iwp) :: file_id_in = 18 !< file unit for input binaray file … … 189 196 190 197 LOGICAL :: convert_average_data = .FALSE. !< namelist parameter to decide whether average or instantaneous data should be converted 198 LOGICAL :: convert_to_netcdf = .FALSE. !< namelist parameter indicating that binary data is also converted to NetCDF files 191 199 LOGICAL, DIMENSION(:), ALLOCATABLE :: eof !< flag indicating that end of binary file is reached 192 200 … … 217 225 !-- VTK format. 218 226 OPEN ( file_id_out_header, FILE = 'HEADER', STATUS = 'REPLACE', & 219 FORM = 'FORMATTED' ) 227 FORM = 'FORMATTED' ) 220 228 ! 221 229 !-- READ grid setup, i.e. the number and position of vertices and surface elements … … 223 231 !-- header information of the VTK file. 224 232 !-- Note, PARAVIEW expects one VTK file for each variable and each timestep. 225 !-- Hence, header information needs to be duplicated multiple times, which will be 226 !-- be done later in a bash script. 233 !-- Hence, header information needs to be duplicated multiple times. 227 234 !-- Moreover, Paraview expects consecutive vertex and polygon data, which are 228 235 !-- all distributed over the binaray files. Hence, first read vertex data from … … 245 252 READ ( file_id_in ) ns(f) 246 253 READ ( file_id_in ) ns_total 247 254 ! 255 !-- If binary surface dater shall be converted into NetCDF, create the file. 256 !-- Only for first call. 257 IF ( convert_to_netcdf .AND. f == 0 ) THEN 258 CALL surface_output_create_netcdf_file 259 ENDIF 248 260 ! 249 261 !-- Allocate arrays where all the surfaces and vertices will be stored. 250 262 ALLOCATE( points(3,1:npoints(f)) ) 251 252 263 ! 253 264 !-- Read polygon data and store them in a temporary file. … … 404 415 ENDIF 405 416 ! 406 !-- Copy HEADER file 407 CALL system('cp HEADER ' // TRIM( path ) // TRIM( char_dum ) // & 408 's_' // TRIM( variable_name ) // '.vtk') 409 !-- Open VTK file. 410 OPEN ( file_id_out, FILE = TRIM( path ) // TRIM( char_dum ) // & 411 's_' // TRIM( variable_name ) // '.vtk', FORM='FORMATTED', & 412 POSITION = 'APPEND' ) 417 !-- Copy HEADER file and open VTK file 418 IF ( convert_average_data ) THEN 419 CALL system('cp HEADER ' // & 420 TRIM( path ) // TRIM( char_dum ) // & 421 '_AV_' // 's_' // TRIM( variable_name ) // '.vtk' ) 422 423 OPEN ( file_id_out, FILE = TRIM( path ) // TRIM( char_dum ) // & 424 '_AV_' // 's_' // TRIM( variable_name ) // '.vtk', & 425 FORM='FORMATTED', POSITION = 'APPEND' ) 426 ELSE 427 CALL system('cp HEADER ' // & 428 TRIM( path ) // TRIM( char_dum ) // & 429 's_' // TRIM( variable_name ) // '.vtk' ) 430 431 OPEN ( file_id_out, FILE = TRIM( path ) // TRIM( char_dum ) // & 432 's_' // TRIM( variable_name ) // '.vtk', & 433 FORM='FORMATTED', POSITION = 'APPEND' ) 434 ENDIF 435 413 436 ENDIF 414 437 ! … … 459 482 INTEGER(iwp) :: file_id_parin = 90 460 483 461 NAMELIST /surface_output/ convert_average_data, c ycle_number, num_pe,&462 path, run484 NAMELIST /surface_output/ convert_average_data, convert_to_netcdf, & 485 cycle_number, num_pe, path, run 463 486 464 487 ! … … 534 557 END SUBROUTINE surface_output_create_file_string 535 558 536 END PROGRAM surface_output_to_vtk 537 538 539 559 560 !------------------------------------------------------------------------------! 561 ! Description: 562 ! ------------ 563 !> This subroutine creates the NetCDF file and defines the dimensions 564 !------------------------------------------------------------------------------! 565 SUBROUTINE surface_output_create_netcdf_file 566 567 IMPLICIT NONE 568 569 570 CHARACTER(LEN=5) :: char_cycle = '' 571 CHARACTER(LEN=200) :: nc_filename 572 573 ! 574 !-- Create substring for the cycle number. 575 IF ( cycle_number /= 0 ) THEN 576 IF ( cycle_number < 10 ) THEN 577 WRITE( char_cycle, '(I1)') cycle_number 578 char_cycle = '.00' // TRIM( char_cycle ) // '.' 579 ELSEIF ( cycle_number < 100 ) THEN 580 WRITE( char_cycle, '(I2)') cycle_number 581 char_cycle = '.0' // TRIM( char_cycle ) // '.' 582 ELSEIF ( cycle_number < 1000 ) THEN 583 WRITE( char_cycle, '(I3)') cycle_number 584 char_cycle = '.' // TRIM( char_cycle ) // '.' 585 ENDIF 586 ELSE 587 char_cycle = '.' 588 ENDIF 589 #if defined( __netcdf ) 590 ! 591 !- Create file 592 nc_filename = TRIM( path ) // TRIM(run) // '_surf' // & 593 TRIM( char_cycle ) // 'nc' 594 595 status_nc = NF90_CREATE( nc_filename, IOR( NF90_CLOBBER, NF90_NETCDF4 ),& 596 nc_id ) 597 CALL handle_error( 'create netcdf file' ) 598 ! 599 !-- Define dimensions 600 status_nc = NF90_DEF_DIM( nc_id, 'xs', ns_total, id_xs ) 601 CALL handle_error( 'define dimension xs' ) 602 status_nc = NF90_DEF_DIM( nc_id, 'ys', ns_total, id_ys ) 603 CALL handle_error( 'define dimension ys' ) 604 status_nc = NF90_DEF_DIM( nc_id, 'zs', ns_total, id_zs ) 605 CALL handle_error( 'define dimension zs' ) 606 status_nc = NF90_DEF_DIM( nc_id, 'azimuth', ns_total, id_az ) 607 CALL handle_error( 'define dimension azimuth' ) 608 status_nc = NF90_DEF_DIM( nc_id, 'zenith', ns_total, id_ze ) 609 CALL handle_error( 'define dimension zenith' ) 610 status_nc = NF90_DEF_DIM( nc_id, 'time', NF90_UNLIMITED, id_time ) 611 CALL handle_error( 'define dimension time' ) 612 ! 613 !-- Define dimension variables 614 status_nc = NF90_DEF_VAR( nc_id, 'xs', NF90_FLOAT, & 615 (/ id_xs /), id_var_xs ) 616 CALL handle_error( 'define variable xs' ) 617 status_nc = NF90_DEF_VAR( nc_id, 'ys', NF90_FLOAT, & 618 (/ id_ys /), id_var_ys ) 619 CALL handle_error( 'define variable ys' ) 620 status_nc = NF90_DEF_VAR( nc_id, 'zs', NF90_FLOAT, & 621 (/ id_zs /), id_var_zs ) 622 CALL handle_error( 'define variable zs' ) 623 status_nc = NF90_DEF_VAR( nc_id, 'azimuth', NF90_FLOAT, & 624 (/ id_az /), id_var_az ) 625 CALL handle_error( 'define variable azimuth' ) 626 status_nc = NF90_DEF_VAR( nc_id, 'zenith', NF90_FLOAT, & 627 (/ id_ze /), id_var_ze ) 628 CALL handle_error( 'define variable zenith' ) 629 status_nc = NF90_DEF_VAR( nc_id, 'time', NF90_FLOAT, & 630 (/ id_time /), id_var_time ) 631 CALL handle_error( 'define variable time' ) 632 #endif 633 END SUBROUTINE surface_output_create_netcdf_file 634 635 636 !------------------------------------------------------------------------------! 637 ! Description: 638 ! ------------ 639 !> Output of NetCDF error code 640 !------------------------------------------------------------------------------! 641 SUBROUTINE handle_error( action ) 642 643 IMPLICIT NONE 644 645 CHARACTER(LEN=*) :: action 646 647 #if defined( __netcdf ) 648 IF ( status_nc /= NF90_NOERR ) THEN 649 PRINT*, TRIM( NF90_STRERROR( status_nc ) ) // ' -- ' // action 650 STOP 651 ENDIF 652 #endif 653 654 END SUBROUTINE handle_error 655 656 END PROGRAM surface_output_to_vtk_netcdf 657 658 659
Note: See TracChangeset
for help on using the changeset viewer.
