!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! ! MODULE OUTPUT_MODULE ! ! This module handles the output of the fields that are generated by the main ! geogrid routines. This output may include output to a console and output to ! the WRF I/O API. !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! module output_module use parallel_module use gridinfo_module use misc_definitions_module use module_debug use dio use wrfheader integer, parameter :: MAX_DIMENSIONS = 3 #ifdef _GEOGRID ! Information about fields that will be written integer :: NUM_AUTOMATIC_FIELDS ! Set later, but very near to a parameter #endif integer :: NUM_FIELDS type field_info integer :: ndims, istagger integer, dimension(MAX_DIMENSIONS) :: dom_start, mem_start, patch_start integer, dimension(MAX_DIMENSIONS) :: dom_end, mem_end, patch_end real, pointer, dimension(:,:,:) :: rdata_arr character (len=128), dimension(MAX_DIMENSIONS) :: dimnames character (len=128) :: fieldname, mem_order, stagger, units, descr end type field_info type (field_info), pointer, dimension(:) :: fields type(dio_file) :: dfile contains !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! ! Name: output_init ! ! Purpose: To initialize the output module. Such initialization may include ! opening an X window, and making initialization calls to an I/O API. !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! subroutine output_init(nest_number, title, datestr, grid_type, dynopt, & corner_lats, corner_lons, & start_dom_1, end_dom_1, start_dom_2, end_dom_2, & start_patch_1, end_patch_1, start_patch_2, end_patch_2, & start_mem_1, end_mem_1, start_mem_2, end_mem_2, & extra_col, extra_row) #ifdef _GEOGRID use llxy_module use source_data_module #endif implicit none ! Arguments integer, intent(in) :: nest_number, dynopt, & start_dom_1, end_dom_1, start_dom_2, end_dom_2, & start_patch_1, end_patch_1, start_patch_2, end_patch_2, & start_mem_1, end_mem_1, start_mem_2, end_mem_2 real, dimension(16), intent(in) :: corner_lats, corner_lons logical, intent(in) :: extra_col, extra_row character (len=1), intent(in) :: grid_type character (len=19), intent(in) :: datestr character (len=*), intent(in) :: title ! Local variables integer :: i, istatus, save_domain, comm_1, comm_2 integer :: sp1, ep1, sp2, ep2, ep1_stag, ep2_stag, N real :: dx, dy, cen_lat, cen_lon, moad_cen_lat, dx_tmp(20), dy_tmp(20) character (len=128) :: coption character (len=MAX_FILENAME_LEN) :: output_fname call init_output_fields(grid_type, & start_dom_1, end_dom_1, start_dom_2, end_dom_2, & start_patch_1, end_patch_1, start_patch_2, end_patch_2, & start_mem_1, end_mem_1, start_mem_2, end_mem_2, & extra_col, extra_row) if (my_proc_id == IO_NODE .or. do_tiled_output) then istatus = 0 call dio_init(iret=istatus) call mprintf((istatus /= 0), ERROR, 'Error in ext_pkg_ioinit') istatus = 0 #ifdef _GEOGRID output_fname = ' ' if (grid_type == 'C') then output_fname = trim(opt_output_from_geogrid_path)//'geo_em.d .dio' i = len_trim(opt_output_from_geogrid_path) write(output_fname(i+9:i+10),'(i2.2)') nest_number else if (grid_type == 'E') then if (nest_number == 1) then output_fname = trim(opt_output_from_geogrid_path)//'geo_nmm.d .dio' i = len_trim(opt_output_from_geogrid_path) write(output_fname(i+10:i+11),'(i2.2)') nest_number else output_fname = trim(opt_output_from_geogrid_path)//'geo_nmm_nest.l .dio' i = len_trim(opt_output_from_geogrid_path) write(output_fname(i+15:i+16),'(i2.2)') nest_number-1 end if else if (grid_type == 'B') then output_fname = trim(opt_output_from_geogrid_path)//'geo_nmb.d .dio' i = len_trim(opt_output_from_geogrid_path) write(output_fname(i+10:i+11),'(i2.2)') nest_number else if (grid_type == 'A') then output_fname = trim(opt_output_from_geogrid_path)//'geo_sla.d .dio' i = len_trim(opt_output_from_geogrid_path) write(output_fname(i+10:i+11),'(i2.2)') nest_number end if if (nprocs > 1 .and. do_tiled_output) then write(output_fname(len_trim(output_fname)+1:len_trim(output_fname)+5), '(a1,i4.4)') & '_', my_proc_id end if #endif #ifdef _METGRID output_fname = ' ' if (grid_type == 'C') then output_fname = trim(opt_output_from_metgrid_path)//'met_em.d .'//trim(datestr)//'.dio' i = len_trim(opt_output_from_metgrid_path) write(output_fname(i+9:i+10),'(i2.2)') nest_number else if (grid_type == 'E') then output_fname = trim(opt_output_from_metgrid_path)//'met_nmm.d .'//trim(datestr)//'.dio' i = len_trim(opt_output_from_metgrid_path) write(output_fname(i+10:i+11),'(i2.2)') nest_number else if (grid_type == 'B') then output_fname = trim(opt_output_from_metgrid_path)//'met_nmb.d .'//trim(datestr)//'.dio' i = len_trim(opt_output_from_metgrid_path) write(output_fname(i+10:i+11),'(i2.2)') nest_number else if (grid_type == 'A') then output_fname = trim(opt_output_from_metgrid_path)//'met_sla.d .'//trim(datestr)//'.dio' i = len_trim(opt_output_from_metgrid_path) write(output_fname(i+10:i+11),'(i2.2)') nest_number end if if (nprocs > 1 .and. do_tiled_output) then write(output_fname(len_trim(output_fname)+1:len_trim(output_fname)+5), '(a1,i4.4)') & '_', my_proc_id end if #endif end if if (my_proc_id == IO_NODE .or. do_tiled_output) then istatus = 0 call dio_open(dfile,trim(output_fname),"WRITE",iret=istatus) call mprintf((istatus /= 0),ERROR,'Error in ext_pkg_open_for_write_begin') end if #ifdef _GEOGRID sp1 = start_patch_1 ep1 = end_patch_1 sp2 = start_patch_2 ep2 = end_patch_2 if (grid_type == 'C') then if (extra_col .or. (my_proc_id == IO_NODE .and. .not. do_tiled_output)) then ep1_stag = ep1 + 1 else ep1_stag = ep1 end if if (extra_row .or. (my_proc_id == IO_NODE .and. .not. do_tiled_output)) then ep2_stag = ep2 + 1 else ep2_stag = ep2 end if else if (grid_type == 'E' .or. grid_type =='B' .or. grid_type == 'A') then ep1 = ep1 ep2 = ep2 ep1_stag = ep1 ep2_stag = ep2 end if if (grid_type == 'C' .or. grid_type == 'A') then dx = proj_stack(nest_number)%dx dy = proj_stack(nest_number)%dx save_domain = iget_selected_domain() ! Note: In the following, we use ixdim/2 rather than (ixdim+1)/2 because ! the i/j coordinates given to xytoll must be with respect to the ! mass grid, and ixdim and jydim are the full grid dimensions. ! Get MOAD_CEN_LAT call select_domain(1) call xytoll(real(ixdim(1))/2.,real(jydim(1))/2., moad_cen_lat, cen_lon, M) ! Get CEN_LAT and CEN_LON for this nest call select_domain(nest_number) call xytoll(real(ixdim(nest_number))/2.,real(jydim(nest_number))/2., cen_lat, cen_lon, M) call select_domain(save_domain) else if (grid_type == 'E') then dx = dxkm / 3**(nest_number-1) ! For NMM, nest_number is really nesting level dy = dykm / 3**(nest_number-1) moad_cen_lat = 0. cen_lat=known_lat cen_lon=known_lon else if (grid_type == 'B') then dx_tmp(1)=dxkm dy_tmp(1)=dykm DO N=2,nest_number dx_tmp(N)=dx_tmp(parent_id(N))/(parent_grid_ratio(N)) dy_tmp(N)=dy_tmp(parent_id(N))/(parent_grid_ratio(N)) ENDDO print*, 'writing dx, dy: ', dx_tmp(nest_number), dy_tmp(nest_number) dx = dx_tmp(nest_number) ! For NMMB, need to worry about parent and ratios dy = dy_tmp(nest_number) moad_cen_lat = 0. cen_lat=known_lat cen_lon=known_lon end if ! We may now write global attributes to the file call write_global_attrs(title, datestr, grid_type, dynopt, ixdim(nest_number), jydim(nest_number), 0, & sp1, ep1, sp1, ep1_stag, sp2, ep2, sp2, ep2_stag, & iproj_type, source_mminlu, source_iswater, source_isice, source_isurban, & source_isoilwater, nest_number, parent_id(nest_number), & proj_stack(nest_number)%i_parent_start,proj_stack(nest_number)%j_parent_start, & proj_stack(nest_number)%i_parent_end ,proj_stack(nest_number)%j_parent_end, & dx, dy, cen_lat, moad_cen_lat, & cen_lon, stand_lon, truelat1, truelat2, & parent_grid_ratio(nest_number), corner_lats, corner_lons) #endif end subroutine output_init !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! ! Name: init_output_fields ! ! Purpose: To fill in structures describing each of the fields that will be ! written to the I/O API !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! subroutine init_output_fields(grid_type, & start_dom_1, end_dom_1, start_dom_2, end_dom_2, & start_patch_1, end_patch_1, start_patch_2, end_patch_2, & start_mem_1, end_mem_1, start_mem_2, end_mem_2, & extra_col, extra_row) ! Modules #ifdef _GEOGRID use source_data_module #endif #ifdef _METGRID use storage_module #endif use parallel_module implicit none ! Arguments integer, intent(in) :: start_dom_1, end_dom_1, start_dom_2, end_dom_2, & start_patch_1, end_patch_1, start_patch_2, end_patch_2, & start_mem_1, end_mem_1, start_mem_2, end_mem_2 logical, intent(in) :: extra_col, extra_row character (len=1), intent(in) :: grid_type ! Local variables integer :: i, istagger, ifieldstatus, & nfields, min_category, max_category integer :: ndims character (len=128) :: fieldname character (len=128) :: memorder, units, description character (len=128), dimension(3) :: dimnames ! ! First find out how many fields there are ! call reset_next_field() ifieldstatus = 0 nfields = 0 do while (ifieldstatus == 0) call get_next_output_fieldname(fieldname, ndims, & min_category, max_category, & istagger, memorder, dimnames, & units, description, ifieldstatus) if (ifieldstatus == 0) nfields = nfields + 1 end do #ifdef _METGRID NUM_FIELDS = nfields #endif #ifdef _GEOGRID if (grid_type == 'C') NUM_AUTOMATIC_FIELDS = 20 if (grid_type == 'E' .or. grid_type == 'B') NUM_AUTOMATIC_FIELDS = 7 if (grid_type == 'A') NUM_AUTOMATIC_FIELDS = 4 NUM_FIELDS = nfields+NUM_AUTOMATIC_FIELDS allocate(fields(NUM_FIELDS)) ! ! There are some fields that will always be computed ! Initialize those fields first, followed by all user-specified fields ! if (grid_type == 'C') then fields(1)%fieldname = 'XLAT_M' fields(1)%units = 'degrees latitude' fields(1)%descr = 'Latitude on mass grid' fields(2)%fieldname = 'XLONG_M' fields(2)%units = 'degrees longitude' fields(2)%descr = 'Longitude on mass grid' fields(3)%fieldname = 'XLAT_V' fields(3)%units = 'degrees latitude' fields(3)%descr = 'Latitude on V grid' fields(4)%fieldname = 'XLONG_V' fields(4)%units = 'degrees longitude' fields(4)%descr = 'Longitude on V grid' fields(5)%fieldname = 'XLAT_U' fields(5)%units = 'degrees latitude' fields(5)%descr = 'Latitude on U grid' fields(6)%fieldname = 'XLONG_U' fields(6)%units = 'degrees longitude' fields(6)%descr = 'Longitude on U grid' fields(7)%fieldname = 'MAPFAC_M' fields(7)%units = 'none' fields(7)%descr = 'Mapfactor on mass grid' fields(8)%fieldname = 'MAPFAC_V' fields(8)%units = 'none' fields(8)%descr = 'Mapfactor on V grid' fields(9)%fieldname = 'MAPFAC_U' fields(9)%units = 'none' fields(9)%descr = 'Mapfactor on U grid' fields(10)%fieldname = 'MAPFAC_MX' fields(10)%units = 'none' fields(10)%descr = 'Mapfactor (x-dir) on mass grid' fields(11)%fieldname = 'MAPFAC_VX' fields(11)%units = 'none' fields(11)%descr = 'Mapfactor (x-dir) on V grid' fields(12)%fieldname = 'MAPFAC_UX' fields(12)%units = 'none' fields(12)%descr = 'Mapfactor (x-dir) on U grid' fields(13)%fieldname = 'MAPFAC_MY' fields(13)%units = 'none' fields(13)%descr = 'Mapfactor (y-dir) on mass grid' fields(14)%fieldname = 'MAPFAC_VY' fields(14)%units = 'none' fields(14)%descr = 'Mapfactor (y-dir) on V grid' fields(15)%fieldname = 'MAPFAC_UY' fields(15)%units = 'none' fields(15)%descr = 'Mapfactor (y-dir) on U grid' fields(16)%fieldname = 'E' fields(16)%units = '-' fields(16)%descr = 'Coriolis E parameter' fields(17)%fieldname = 'F' fields(17)%units = '-' fields(17)%descr = 'Coriolis F parameter' fields(18)%fieldname = 'SINALPHA' fields(18)%units = 'none' fields(18)%descr = 'Sine of rotation angle' fields(19)%fieldname = 'COSALPHA' fields(19)%units = 'none' fields(19)%descr = 'Cosine of rotation angle' fields(20)%fieldname = 'LANDMASK' fields(20)%units = 'none' fields(20)%descr = 'Landmask : 1=land, 0=water' else if (grid_type == 'E' .or. grid_type == 'B') then fields(1)%fieldname = 'XLAT_M' fields(1)%units = 'degrees latitude' fields(1)%descr = 'Latitude on mass grid' fields(2)%fieldname = 'XLONG_M' fields(2)%units = 'degrees longitude' fields(2)%descr = 'Longitude on mass grid' fields(3)%fieldname = 'XLAT_V' fields(3)%units = 'degrees latitude' fields(3)%descr = 'Latitude on velocity grid' fields(4)%fieldname = 'XLONG_V' fields(4)%units = 'degrees longitude' fields(4)%descr = 'Longitude on velocity grid' fields(5)%fieldname = 'E' fields(5)%units = '-' fields(5)%descr = 'Coriolis E parameter' fields(6)%fieldname = 'F' fields(6)%units = '-' fields(6)%descr = 'Coriolis F parameter' fields(7)%fieldname = 'LANDMASK' fields(7)%units = 'none' fields(7)%descr = 'Landmask : 1=land, 0=water' else if (grid_type == 'A') then fields(1)%fieldname = 'XLAT_M' fields(1)%units = 'degrees latitude' fields(1)%descr = 'Latitude on mass grid' fields(2)%fieldname = 'XLONG_M' fields(2)%units = 'degrees longitude' fields(2)%descr = 'Longitude on mass grid' fields(3)%fieldname = 'MAPFAC_M' fields(3)%units = 'none' fields(3)%descr = 'Mapfactor on mass grid' ! fields(4)%fieldname = 'E' ! fields(4)%units = '-' ! fields(4)%descr = 'Coriolis E parameter' ! fields(5)%fieldname = 'F' ! fields(5)%units = '-' ! fields(5)%descr = 'Coriolis F parameter' fields(4)%fieldname = 'LANDMASK' fields(4)%units = 'none' fields(4)%descr = 'Landmask : 1=land, 0=water' end if ! ! General defaults for "always computed" fields ! do i=1,NUM_AUTOMATIC_FIELDS fields(i)%ndims = 2 fields(i)%dom_start(1) = start_dom_1 fields(i)%dom_start(2) = start_dom_2 fields(i)%dom_start(3) = 1 fields(i)%mem_start(1) = start_mem_1 fields(i)%mem_start(2) = start_mem_2 fields(i)%mem_start(3) = 1 fields(i)%patch_start(1) = start_patch_1 fields(i)%patch_start(2) = start_patch_2 fields(i)%patch_start(3) = 1 fields(i)%dom_end(1) = end_dom_1 fields(i)%dom_end(2) = end_dom_2 fields(i)%dom_end(3) = 1 fields(i)%mem_end(1) = end_mem_1 fields(i)%mem_end(2) = end_mem_2 fields(i)%mem_end(3) = 1 fields(i)%patch_end(1) = end_patch_1 fields(i)%patch_end(2) = end_patch_2 fields(i)%patch_end(3) = 1 fields(i)%dimnames(3) = ' ' fields(i)%mem_order = 'XY' fields(i)%stagger = 'M' if (grid_type == 'C' .or. grid_type == 'A') then fields(i)%istagger = M else if (grid_type == 'E' .or. grid_type == 'B') then fields(i)%istagger = HH end if fields(i)%dimnames(1) = 'west_east' fields(i)%dimnames(2) = 'south_north' end do if (grid_type == 'C') then ! Lat V if (extra_row) then fields(3)%dom_end(2) = fields(3)%dom_end(2) + 1 fields(3)%mem_end(2) = fields(3)%mem_end(2) + 1 fields(3)%patch_end(2) = fields(3)%patch_end(2) + 1 else if (my_proc_id == IO_NODE .and. .not. do_tiled_output) then fields(3)%dom_end(2) = fields(3)%dom_end(2) + 1 end if fields(3)%dimnames(2) = 'south_north_stag' fields(3)%stagger = 'V' fields(3)%istagger = V ! Lon V if (extra_row) then fields(4)%dom_end(2) = fields(4)%dom_end(2) + 1 fields(4)%mem_end(2) = fields(4)%mem_end(2) + 1 fields(4)%patch_end(2) = fields(4)%patch_end(2) + 1 else if (my_proc_id == IO_NODE .and. .not. do_tiled_output) then fields(4)%dom_end(2) = fields(4)%dom_end(2) + 1 end if fields(4)%dimnames(2) = 'south_north_stag' fields(4)%stagger = 'V' fields(4)%istagger = V ! Lat U if (extra_col) then fields(5)%dom_end(1) = fields(5)%dom_end(1) + 1 fields(5)%mem_end(1) = fields(5)%mem_end(1) + 1 fields(5)%patch_end(1) = fields(5)%patch_end(1) + 1 else if (my_proc_id == IO_NODE .and. .not. do_tiled_output) then fields(5)%dom_end(1) = fields(5)%dom_end(1) + 1 end if fields(5)%dimnames(1) = 'west_east_stag' fields(5)%stagger = 'U' fields(5)%istagger = U ! Lon U if (extra_col) then fields(6)%dom_end(1) = fields(6)%dom_end(1) + 1 fields(6)%mem_end(1) = fields(6)%mem_end(1) + 1 fields(6)%patch_end(1) = fields(6)%patch_end(1) + 1 else if (my_proc_id == IO_NODE .and. .not. do_tiled_output) then fields(6)%dom_end(1) = fields(6)%dom_end(1) + 1 end if fields(6)%dimnames(1) = 'west_east_stag' fields(6)%stagger = 'U' fields(6)%istagger = U ! Mapfac V if (extra_row) then fields(8)%dom_end(2) = fields(8)%dom_end(2) + 1 fields(8)%mem_end(2) = fields(8)%mem_end(2) + 1 fields(8)%patch_end(2) = fields(8)%patch_end(2) + 1 else if (my_proc_id == IO_NODE .and. .not. do_tiled_output) then fields(8)%dom_end(2) = fields(8)%dom_end(2) + 1 end if fields(8)%dimnames(2) = 'south_north_stag' fields(8)%stagger = 'V' fields(8)%istagger = V ! Mapfac U if (extra_col) then fields(9)%dom_end(1) = fields(9)%dom_end(1) + 1 fields(9)%mem_end(1) = fields(9)%mem_end(1) + 1 fields(9)%patch_end(1) = fields(9)%patch_end(1) + 1 else if (my_proc_id == IO_NODE .and. .not. do_tiled_output) then fields(9)%dom_end(1) = fields(9)%dom_end(1) + 1 end if fields(9)%dimnames(1) = 'west_east_stag' fields(9)%stagger = 'U' fields(9)%istagger = U ! Mapfac V-X if (extra_row) then fields(11)%dom_end(2) = fields(11)%dom_end(2) + 1 fields(11)%mem_end(2) = fields(11)%mem_end(2) + 1 fields(11)%patch_end(2) = fields(11)%patch_end(2) + 1 else if (my_proc_id == IO_NODE .and. .not. do_tiled_output) then fields(11)%dom_end(2) = fields(11)%dom_end(2) + 1 end if fields(11)%dimnames(2) = 'south_north_stag' fields(11)%stagger = 'V' fields(11)%istagger = V ! Mapfac U-X if (extra_col) then fields(12)%dom_end(1) = fields(12)%dom_end(1) + 1 fields(12)%mem_end(1) = fields(12)%mem_end(1) + 1 fields(12)%patch_end(1) = fields(12)%patch_end(1) + 1 else if (my_proc_id == IO_NODE .and. .not. do_tiled_output) then fields(12)%dom_end(1) = fields(12)%dom_end(1) + 1 end if fields(12)%dimnames(1) = 'west_east_stag' fields(12)%stagger = 'U' fields(12)%istagger = U ! Mapfac V-Y if (extra_row) then fields(14)%dom_end(2) = fields(14)%dom_end(2) + 1 fields(14)%mem_end(2) = fields(14)%mem_end(2) + 1 fields(14)%patch_end(2) = fields(14)%patch_end(2) + 1 else if (my_proc_id == IO_NODE .and. .not. do_tiled_output) then fields(14)%dom_end(2) = fields(14)%dom_end(2) + 1 end if fields(14)%dimnames(2) = 'south_north_stag' fields(14)%stagger = 'V' fields(14)%istagger = V ! Mapfac U-Y if (extra_col) then fields(15)%dom_end(1) = fields(15)%dom_end(1) + 1 fields(15)%mem_end(1) = fields(15)%mem_end(1) + 1 fields(15)%patch_end(1) = fields(15)%patch_end(1) + 1 else if (my_proc_id == IO_NODE .and. .not. do_tiled_output) then fields(15)%dom_end(1) = fields(15)%dom_end(1) + 1 end if fields(15)%dimnames(1) = 'west_east_stag' fields(15)%stagger = 'U' fields(15)%istagger = U else if (grid_type == 'E' .or. grid_type == 'B') then ! Lat V fields(3)%stagger = 'V' fields(3)%istagger = VV ! Lon V fields(4)%stagger = 'V' fields(4)%istagger = VV end if #endif ! ! Now set up the field_info structure for each user-specified field ! call reset_next_field() ifieldstatus = 0 #ifdef _GEOGRID nfields = NUM_AUTOMATIC_FIELDS+1 #endif #ifdef _METGRID allocate(fields(NUM_FIELDS)) nfields = 1 #endif do while (ifieldstatus == 0) !{ call get_next_output_fieldname(fieldname, ndims, & min_category, max_category, & istagger, memorder, dimnames, & units, description, ifieldstatus) if (ifieldstatus == 0) then !{ fields(nfields)%ndims = ndims fields(nfields)%fieldname = fieldname fields(nfields)%istagger = istagger if (istagger == M) then fields(nfields)%stagger = 'M' else if (istagger == U) then fields(nfields)%stagger = 'U' else if (istagger == V) then fields(nfields)%stagger = 'V' else if (istagger == HH) then fields(nfields)%stagger = 'M' else if (istagger == VV) then fields(nfields)%stagger = 'V' end if fields(nfields)%mem_order = memorder fields(nfields)%dimnames(1) = dimnames(1) fields(nfields)%dimnames(2) = dimnames(2) fields(nfields)%dimnames(3) = dimnames(3) fields(nfields)%units = units fields(nfields)%descr = description fields(nfields)%dom_start(1)=start_dom_1 fields(nfields)%dom_start(2)=start_dom_2 fields(nfields)%dom_start(3)=min_category fields(nfields)%mem_start(1)=start_mem_1 fields(nfields)%mem_start(2)=start_mem_2 fields(nfields)%mem_start(3)=min_category fields(nfields)%patch_start(1)=start_patch_1 fields(nfields)%patch_start(2)=start_patch_2 fields(nfields)%patch_start(3)=min_category fields(nfields)%dom_end(1)=end_dom_1 fields(nfields)%dom_end(2)=end_dom_2 fields(nfields)%dom_end(3)=max_category fields(nfields)%mem_end(1)=end_mem_1 fields(nfields)%mem_end(2)=end_mem_2 fields(nfields)%mem_end(3)=max_category fields(nfields)%patch_end(1)=end_patch_1 fields(nfields)%patch_end(2)=end_patch_2 fields(nfields)%patch_end(3)=max_category if (extra_col .and. istagger == U) then !{ fields(nfields)%dom_end(1)=fields(nfields)%dom_end(1) + 1 fields(nfields)%mem_end(1)=fields(nfields)%mem_end(1) + 1 fields(nfields)%patch_end(1)=fields(nfields)%patch_end(1) + 1 else if (istagger == U .and. my_proc_id == IO_NODE .and. .not. do_tiled_output) then fields(nfields)%dom_end(1)=fields(nfields)%dom_end(1) + 1 end if !} if (extra_row .and. istagger == V) then !{ fields(nfields)%dom_end(2)=fields(nfields)%dom_end(2) + 1 fields(nfields)%mem_end(2)=fields(nfields)%mem_end(2) + 1 fields(nfields)%patch_end(2)=fields(nfields)%patch_end(2) + 1 else if (istagger == V .and. my_proc_id == IO_NODE .and. .not. do_tiled_output) then fields(nfields)%dom_end(2)=fields(nfields)%dom_end(2) + 1 end if !} nfields = nfields + 1 end if ! the next field given by get_next_fieldname() is valid } end do ! for each user-specified field } end subroutine init_output_fields !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! ! Name: write_field ! ! Purpose: This routine writes the provided field to any output devices or APIs !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! subroutine write_field(start_mem_i, end_mem_i, & start_mem_j, end_mem_j, & start_mem_k, end_mem_k, & cname, datestr, real_array, is_training) implicit none ! Arguments integer, intent(in) :: start_mem_i, end_mem_i, start_mem_j, end_mem_j, start_mem_k, end_mem_k real, target, dimension(start_mem_i:end_mem_i, start_mem_j:end_mem_j, start_mem_k:end_mem_k), & intent(in) :: real_array logical, intent(in), optional :: is_training character (len=19), intent(in) :: datestr character (len=*), intent(in) :: cname ! Local variables integer :: i integer :: istatus, comm_1, comm_2, domain_desc integer, dimension(3) :: sd, ed, sp, ep, sm, em real, pointer, dimension(:,:,:) :: real_dom_array logical :: allocated_real_locally integer, dimension(512) :: header allocated_real_locally = .false. ! If we are running distributed memory and need to gather all tiles onto a single processor for output if (nprocs > 1 .and. .not. do_tiled_output) then do i=1,NUM_FIELDS if ( (index(cname, trim(fields(i)%fieldname) ) /= 0) .and. & (len_trim(cname) == len_trim(fields(i)%fieldname)) ) then istatus = 0 ! For the gather routines below, the IO_NODE should give the full domain dimensions, but the ! memory and patch dimensions should indicate what the processor already has in its patch_array. ! This is because an array with dimensions of the full domain will be allocated, and the patch_array ! will be copied from local memory into the full domain array in the area specified by the patch ! dimensions. sd = fields(i)%dom_start ed = fields(i)%dom_end sp = fields(i)%patch_start ep = fields(i)%patch_end sm = fields(i)%mem_start em = fields(i)%mem_end allocate(real_dom_array(sd(1):ed(1),sd(2):ed(2),sd(3):ed(3))) allocated_real_locally = .true. call gather_whole_field_r(real_array, & sm(1), em(1), sm(2), em(2), sm(3), em(3), & sp(1), ep(1), sp(2), ep(2), sp(3), ep(3), & real_dom_array, & sd(1), ed(1), sd(2), ed(2), sd(3), ed(3)) exit end if end do else do i=1,NUM_FIELDS if ( (index(cname, trim(fields(i)%fieldname) ) /= 0) .and. & (len_trim(cname) == len_trim(fields(i)%fieldname)) ) then istatus = 0 real_dom_array => real_array exit end if end do end if ! Now a write call is only done if each processor writes its own file, or if we are the IO_NODE if (my_proc_id == IO_NODE .or. do_tiled_output) then comm_1 = 1 comm_2 = 1 domain_desc = 0 do i=1,NUM_FIELDS if ( (index(cname, trim(fields(i)%fieldname) ) /= 0) .and. & (len_trim(cname) == len_trim(fields(i)%fieldname)) ) then ! Here, the output array has dimensions of the full grid if it was gathered together ! from all processors if (my_proc_id == IO_NODE .and. nprocs > 1 .and. .not. do_tiled_output) then sd = fields(i)%dom_start ed = fields(i)%dom_end sm = sd em = ed sp = sd ep = ed ! If we are writing one file per processor, then each processor only writes out the ! part of the domain that it has in memory else ! BUG: Shouldn't we set sd/ed to be domain_start/domain_end? ! Maybe not, since patch is already adjusted for staggering; but maybe so, and also adjust ! for staggering if it is alright to pass true domain dimensions to write_field. sd = fields(i)%patch_start ed = fields(i)%patch_end sp = fields(i)%patch_start ep = fields(i)%patch_end sm = fields(i)%mem_start em = fields(i)%mem_end end if istatus = 0 header=0 call wrfheader_pack(header, DateStr=datestr, & VarName=trim(fields(i)%fieldname), & Units=trim(fields(i)%units), & Description=trim(fields(i)%descr), & FieldType=DIO_REAL4, & MemoryOrder=trim(fields(i)%mem_order), & Stagger=trim(fields(i)%stagger), & DimNames=fields(i)%dimnames, & DomainStart=sd, DomainEnd=ed, & PatchStart=sp, PatchEnd=ep & ) call dio_write(dfile,trim(fields(i)%fieldname), & real_dom_array(sd(1):ed(1),sd(2):ed(2),sd(3):ed(3)), & header=header, iret=istatus) call mprintf((istatus /= 0),ERROR,'Error in ext_pkg_write_field') exit end if end do end if if (allocated_real_locally) deallocate(real_dom_array) end subroutine write_field !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! ! Name: write_global_attrs ! ! Purpose: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! subroutine write_global_attrs(title, start_date, grid_type, dyn_opt, & west_east_dim, south_north_dim, bottom_top_dim, & we_patch_s, we_patch_e, we_patch_s_stag, we_patch_e_stag, & sn_patch_s, sn_patch_e, sn_patch_s_stag, sn_patch_e_stag, & map_proj, cmminlu, is_water, is_ice, is_urban, i_soilwater, & grid_id, parent_id, i_parent_start, j_parent_start, & i_parent_end, j_parent_end, dx, dy, cen_lat, moad_cen_lat, cen_lon, & stand_lon, truelat1, truelat2, parent_grid_ratio, corner_lats, corner_lons, & flags, nflags) implicit none ! Arguments integer, intent(in) :: dyn_opt, west_east_dim, south_north_dim, bottom_top_dim, & we_patch_s, we_patch_e, we_patch_s_stag, we_patch_e_stag, & sn_patch_s, sn_patch_e, sn_patch_s_stag, sn_patch_e_stag, & map_proj, is_water, is_ice, is_urban, i_soilwater, & grid_id, parent_id, i_parent_start, j_parent_start, & i_parent_end, j_parent_end, parent_grid_ratio integer, intent(in), optional :: nflags real, intent(in) :: dx, dy, cen_lat, moad_cen_lat, cen_lon, stand_lon, truelat1, truelat2 real, dimension(16), intent(in) :: corner_lats, corner_lons character (len=*), intent(in) :: title, start_date, grid_type character (len=128), intent(in) :: cmminlu character (len=128), dimension(:), intent(in), optional :: flags ! Local variables integer :: local_we_patch_s, local_we_patch_s_stag, & local_we_patch_e, local_we_patch_e_stag, & local_sn_patch_s, local_sn_patch_s_stag, & local_sn_patch_e, local_sn_patch_e_stag integer :: i real, dimension(16) :: local_corner_lats, local_corner_lons local_we_patch_s = we_patch_s local_we_patch_s_stag = we_patch_s_stag local_we_patch_e = we_patch_e local_we_patch_e_stag = we_patch_e_stag local_sn_patch_s = sn_patch_s local_sn_patch_s_stag = sn_patch_s_stag local_sn_patch_e = sn_patch_e local_sn_patch_e_stag = sn_patch_e_stag local_corner_lats = corner_lats local_corner_lons = corner_lons if (nprocs > 1) then if (.not. do_tiled_output) then call parallel_bcast_int(local_we_patch_s, processors(0, 0)) call parallel_bcast_int(local_we_patch_s_stag, processors(0, 0)) call parallel_bcast_int(local_sn_patch_s, processors(0, 0)) call parallel_bcast_int(local_sn_patch_s_stag, processors(0, 0)) call parallel_bcast_int(local_we_patch_e, processors(nproc_x-1, nproc_y-1)) call parallel_bcast_int(local_we_patch_e_stag, processors(nproc_x-1, nproc_y-1)) call parallel_bcast_int(local_sn_patch_e, processors(nproc_x-1, nproc_y-1)) call parallel_bcast_int(local_sn_patch_e_stag, processors(nproc_x-1, nproc_y-1)) end if call parallel_bcast_real(local_corner_lats(1), processors(0, 0)) call parallel_bcast_real(local_corner_lats(2), processors(0, nproc_y-1)) call parallel_bcast_real(local_corner_lats(3), processors(nproc_x-1, nproc_y-1)) call parallel_bcast_real(local_corner_lats(4), processors(nproc_x-1, 0)) call parallel_bcast_real(local_corner_lats(5), processors(0, 0)) call parallel_bcast_real(local_corner_lats(6), processors(0, nproc_y-1)) call parallel_bcast_real(local_corner_lats(7), processors(nproc_x-1, nproc_y-1)) call parallel_bcast_real(local_corner_lats(8), processors(nproc_x-1, 0)) call parallel_bcast_real(local_corner_lats(9), processors(0, 0)) call parallel_bcast_real(local_corner_lats(10), processors(0, nproc_y-1)) call parallel_bcast_real(local_corner_lats(11), processors(nproc_x-1, nproc_y-1)) call parallel_bcast_real(local_corner_lats(12), processors(nproc_x-1, 0)) call parallel_bcast_real(local_corner_lats(13), processors(0, 0)) call parallel_bcast_real(local_corner_lats(14), processors(0, nproc_y-1)) call parallel_bcast_real(local_corner_lats(15), processors(nproc_x-1, nproc_y-1)) call parallel_bcast_real(local_corner_lats(16), processors(nproc_x-1, 0)) call parallel_bcast_real(local_corner_lons(1), processors(0, 0)) call parallel_bcast_real(local_corner_lons(2), processors(0, nproc_y-1)) call parallel_bcast_real(local_corner_lons(3), processors(nproc_x-1, nproc_y-1)) call parallel_bcast_real(local_corner_lons(4), processors(nproc_x-1, 0)) call parallel_bcast_real(local_corner_lons(5), processors(0, 0)) call parallel_bcast_real(local_corner_lons(6), processors(0, nproc_y-1)) call parallel_bcast_real(local_corner_lons(7), processors(nproc_x-1, nproc_y-1)) call parallel_bcast_real(local_corner_lons(8), processors(nproc_x-1, 0)) call parallel_bcast_real(local_corner_lons(9), processors(0, 0)) call parallel_bcast_real(local_corner_lons(10), processors(0, nproc_y-1)) call parallel_bcast_real(local_corner_lons(11), processors(nproc_x-1, nproc_y-1)) call parallel_bcast_real(local_corner_lons(12), processors(nproc_x-1, 0)) call parallel_bcast_real(local_corner_lons(13), processors(0, 0)) call parallel_bcast_real(local_corner_lons(14), processors(0, nproc_y-1)) call parallel_bcast_real(local_corner_lons(15), processors(nproc_x-1, nproc_y-1)) call parallel_bcast_real(local_corner_lons(16), processors(nproc_x-1, 0)) end if if (my_proc_id == IO_NODE .or. do_tiled_output) then call ext_put_dom_ti_char ('TITLE', title) call ext_put_dom_ti_char ('SIMULATION_START_DATE', start_date) call ext_put_dom_ti_integer_scalar('WEST-EAST_GRID_DIMENSION', west_east_dim) call ext_put_dom_ti_integer_scalar('SOUTH-NORTH_GRID_DIMENSION', south_north_dim) call ext_put_dom_ti_integer_scalar('BOTTOM-TOP_GRID_DIMENSION', bottom_top_dim) call ext_put_dom_ti_integer_scalar('WEST-EAST_PATCH_START_UNSTAG', local_we_patch_s) call ext_put_dom_ti_integer_scalar('WEST-EAST_PATCH_END_UNSTAG', local_we_patch_e) call ext_put_dom_ti_integer_scalar('WEST-EAST_PATCH_START_STAG', local_we_patch_s_stag) call ext_put_dom_ti_integer_scalar('WEST-EAST_PATCH_END_STAG', local_we_patch_e_stag) call ext_put_dom_ti_integer_scalar('SOUTH-NORTH_PATCH_START_UNSTAG', local_sn_patch_s) call ext_put_dom_ti_integer_scalar('SOUTH-NORTH_PATCH_END_UNSTAG', local_sn_patch_e) call ext_put_dom_ti_integer_scalar('SOUTH-NORTH_PATCH_START_STAG', local_sn_patch_s_stag) call ext_put_dom_ti_integer_scalar('SOUTH-NORTH_PATCH_END_STAG', local_sn_patch_e_stag) call ext_put_dom_ti_char ('GRIDTYPE', grid_type) call ext_put_dom_ti_real_scalar ('DX', dx) call ext_put_dom_ti_real_scalar ('DY', dy) call ext_put_dom_ti_integer_scalar('DYN_OPT', dyn_opt) call ext_put_dom_ti_real_scalar ('CEN_LAT', cen_lat) call ext_put_dom_ti_real_scalar ('CEN_LON', cen_lon) call ext_put_dom_ti_real_scalar ('TRUELAT1', truelat1) call ext_put_dom_ti_real_scalar ('TRUELAT2', truelat2) call ext_put_dom_ti_real_scalar ('MOAD_CEN_LAT', moad_cen_lat) call ext_put_dom_ti_real_scalar ('STAND_LON', stand_lon) call ext_put_dom_ti_real_vector ('corner_lats', local_corner_lats, 16) call ext_put_dom_ti_real_vector ('corner_lons', local_corner_lons, 16) call ext_put_dom_ti_integer_scalar('MAP_PROJ', map_proj) call ext_put_dom_ti_char ('MMINLU', trim(cmminlu)) call ext_put_dom_ti_integer_scalar('ISWATER', is_water) call ext_put_dom_ti_integer_scalar('ISICE', is_ice) call ext_put_dom_ti_integer_scalar('ISURBAN', is_urban) call ext_put_dom_ti_integer_scalar('ISOILWATER', i_soilwater) call ext_put_dom_ti_integer_scalar('grid_id', grid_id) call ext_put_dom_ti_integer_scalar('parent_id', parent_id) call ext_put_dom_ti_integer_scalar('i_parent_start', i_parent_start) call ext_put_dom_ti_integer_scalar('j_parent_start', j_parent_start) call ext_put_dom_ti_integer_scalar('i_parent_end', i_parent_end) call ext_put_dom_ti_integer_scalar('j_parent_end', j_parent_end) call ext_put_dom_ti_integer_scalar('parent_grid_ratio', parent_grid_ratio) #ifdef _METGRID call ext_put_dom_ti_integer_scalar('FLAG_METGRID', 1) #endif if (present(nflags) .and. present(flags)) then do i=1,nflags if (flags(i) /= ' ') then call ext_put_dom_ti_integer_scalar(trim(flags(i)), 1) end if end do end if end if end subroutine write_global_attrs !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! ! Name: ext_put_dom_ti_integer ! ! Purpose: Write a domain time-independent integer attribute to output. !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! subroutine ext_put_dom_ti_integer_scalar(var_name, var_value) implicit none ! Arguments integer, intent(in) :: var_value character (len=*), intent(in) :: var_name ! Local variables integer :: istatus call dio_write(dfile,trim(var_name),var_value,iret=istatus) call mprintf((istatus /= 0),ERROR,'Error in writing domain time-independent attribute') end subroutine ext_put_dom_ti_integer_scalar !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! ! Name: ext_put_dom_ti_integer ! ! Purpose: Write a domain time-independent integer attribute to output. !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! subroutine ext_put_dom_ti_integer_vector(var_name, var_value, n) implicit none ! Arguments integer, intent(in) :: n integer, dimension(n), intent(in) :: var_value character (len=*), intent(in) :: var_name ! Local variables integer :: istatus call dio_write(dfile,trim(var_name),var_value(1:n),iret=istatus) call mprintf((istatus /= 0),ERROR,'Error in writing domain time-independent attribute') end subroutine ext_put_dom_ti_integer_vector !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! ! Name: ext_put_dom_ti_real ! ! Purpose: Write a domain time-independent real attribute to output. !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! subroutine ext_put_dom_ti_real_scalar(var_name, var_value) implicit none ! Arguments real, intent(in) :: var_value character (len=*), intent(in) :: var_name ! Local variables integer :: istatus call dio_write(dfile,trim(var_name),var_value,iret=istatus) call mprintf((istatus /= 0),ERROR,'Error in writing domain time-independent attribute') end subroutine ext_put_dom_ti_real_scalar !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! ! Name: ext_put_dom_ti_real ! ! Purpose: Write a domain time-independent real attribute to output. !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! subroutine ext_put_dom_ti_real_vector(var_name, var_value, n) implicit none ! Arguments integer, intent(in) :: n real, dimension(n), intent(in) :: var_value character (len=*), intent(in) :: var_name ! Local variables integer :: istatus call dio_write(dfile,trim(var_name),var_value(1:n),iret=istatus) call mprintf((istatus /= 0),ERROR,'Error in writing domain time-independent attribute') end subroutine ext_put_dom_ti_real_vector !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! ! Name: ext_put_dom_ti_char ! ! Purpose: Write a domain time-independent character attribute to output. !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! subroutine ext_put_dom_ti_char(var_name, var_value) implicit none ! Arguments character (len=*), intent(in) :: var_name, var_value ! Local variables integer :: istatus call dio_write(dfile,trim(var_name),var_value,iret=istatus) call mprintf((istatus /= 0),ERROR,'Error in writing domain time-independent attribute') end subroutine ext_put_dom_ti_char !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! ! Name: output_close ! ! Purpose: Finalizes all output. This may include closing windows, calling I/O ! API termination routines, or closing files. !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! subroutine output_close() implicit none ! Local variables integer :: istatus if (my_proc_id == IO_NODE .or. do_tiled_output) then istatus = 0 call dio_close(dfile,iret=istatus) call mprintf((istatus /= 0), ERROR, 'Error in ext_pkg_ioclose') istatus = 0 call dio_finalize() call mprintf((istatus /= 0), ERROR, 'Error in ext_pkg_ioexit') end if if (associated(fields)) deallocate(fields) end subroutine output_close end module output_module