subroutine da_write_increments (grid, q_cgrid, mu_cgrid, ph_cgrid) !---------------------------------------------------------------------- ! Purpose: Write analysis increments !---------------------------------------------------------------------- implicit none type (domain), intent(inout) :: grid real,intent(in) :: q_cgrid(ims:ime,jms:jme,kms:kme) real,intent(in) :: ph_cgrid(ims:ime,jms:jme,kms:kme) real,intent(in) :: mu_cgrid(ims:ime,jms:jme) ! Arrays for write out increments: integer :: ix, jy, kz #ifdef DM_PARALLEL real, dimension(1:grid%xb%mix,1:grid%xb%mjy) :: gbuf_2d real, dimension(1:grid%xb%mix+1,1:grid%xb%mjy+1) :: gbuf_2dd real, dimension(1:grid%xb%mix,1:grid%xb%mjy,1:grid%xb%mkz) :: gbuf real, dimension(1:grid%xb%mix,1:grid%xb%mjy,1:grid%xb%mkz+1) :: wgbuf real, dimension(:,:,:), allocatable :: u_global, v_global, w_global, & p_global, t_global, q_global, & ph_global real, dimension(:,:) , allocatable :: mu_global, psfc_global, & psac_global, tgrn_global, terr_global, snow_global,& lat_global, lon_global, lanu_global, & map_factor_global, cori_global, landmask_global #endif integer :: anl_inc_unit if (trace_use) call da_trace_entry("da_write_increments") ! Dimension of the domain: ix = grid%xb%mix jy = grid%xb%mjy kz = grid%xb%mkz #ifdef DM_PARALLEL ! 3-d and 2-d increments: allocate ( p_global (1:ix+1,1:jy+1,1:kz+1)) allocate ( t_global (1:ix+1,1:jy+1,1:kz+1)) allocate ( q_global (1:ix+1,1:jy+1,1:kz+1)) allocate ( u_global (1:ix+1,1:jy+1,1:kz+1)) allocate ( v_global (1:ix+1,1:jy+1,1:kz+1)) allocate ( w_global (1:ix+1,1:jy+1,1:kz+1)) allocate ( ph_global (1:ix+1,1:jy+1,1:kz+1)) allocate (psfc_global (1:ix+1,1:jy+1)) allocate ( mu_global (1:ix+1,1:jy+1)) call da_patch_to_global(grid, grid%xa % p, gbuf) if (rootproc) then p_global(1:ix,1:jy,1:kz) = gbuf(1:ix,1:jy,1:kz) end if call da_patch_to_global(grid, grid%xa % t, gbuf) if (rootproc) then t_global(1:ix,1:jy,1:kz) = gbuf(1:ix,1:jy,1:kz) end if call da_patch_to_global(grid, q_cgrid, gbuf) if (rootproc) then q_global(1:ix,1:jy,1:kz) = gbuf(1:ix,1:jy,1:kz) end if call da_patch_to_global(grid, grid%xa % u, gbuf) if (rootproc) then u_global(1:ix,1:jy,1:kz) = gbuf(1:ix,1:jy,1:kz) end if call da_patch_to_global(grid, grid%xa % v, gbuf) if (rootproc) then v_global(1:ix,1:jy,1:kz) = gbuf(1:ix,1:jy,1:kz) end if ! One more level for w and ph: grid%xp%kde=grid%xp%kde+1 kde=kde+1 call da_patch_to_global(grid, grid%xa % w, wgbuf) if (rootproc) then w_global(1:ix,1:jy,1:kz+1) = wgbuf(1:ix,1:jy,1:kz+1) end if call da_patch_to_global(grid, ph_cgrid, wgbuf) if (rootproc) then ph_global(1:ix,1:jy,1:kz+1) = wgbuf(1:ix,1:jy,1:kz+1) end if kde=kde-1 grid%xp%kde=grid%xp%kde-1 call da_patch_to_global(grid, grid%xa % psfc, gbuf_2d) if (rootproc) then psfc_global(1:ix,1:jy) = gbuf_2d(1:ix,1:jy) end if call da_patch_to_global(grid, mu_cgrid, gbuf_2d) if (rootproc) then mu_global(1:ix,1:jy) = gbuf_2d(1:ix,1:jy) end if ! 2d constant fields: allocate ( psac_global (1:ix+1,1:jy+1)) allocate ( tgrn_global (1:ix+1,1:jy+1)) allocate ( terr_global (1:ix+1,1:jy+1)) allocate ( snow_global (1:ix+1,1:jy+1)) allocate ( lat_global (1:ix+1,1:jy+1)) allocate ( lon_global (1:ix+1,1:jy+1)) allocate ( lanu_global (1:ix+1,1:jy+1)) allocate (map_factor_global (1:ix+1,1:jy+1)) allocate ( cori_global (1:ix+1,1:jy+1)) allocate ( landmask_global (1:ix+1,1:jy+1)) call da_patch_to_global(grid, grid%xb%psac, gbuf_2d) if (rootproc) then psac_global(1:ix,1:jy) = gbuf_2d(1:ix,1:jy) end if call da_patch_to_global(grid, grid%xb%tgrn, gbuf_2d) if (rootproc) then tgrn_global(1:ix,1:jy) = gbuf_2d(1:ix,1:jy) end if call da_patch_to_global(grid, grid%xb%terr, gbuf_2d) if (rootproc) then terr_global(1:ix,1:jy) = gbuf_2d(1:ix,1:jy) end if call da_patch_to_global(grid, grid%xb%snow, gbuf_2d) if (rootproc) then snow_global(1:ix,1:jy) = gbuf_2d(1:ix,1:jy) end if call da_patch_to_global(grid, grid%xb%lat , gbuf_2d) if (rootproc) then lat_global(1:ix,1:jy) = gbuf_2d(1:ix,1:jy) end if call da_patch_to_global(grid, grid%xb%lon , gbuf_2d) if (rootproc) then lon_global(1:ix,1:jy) = gbuf_2d(1:ix,1:jy) end if call da_patch_to_global(grid, grid%xb%lanu, gbuf_2d) if (rootproc) then lanu_global(1:ix,1:jy) = gbuf_2d(1:ix,1:jy) end if call da_patch_to_global(grid, grid%xb%map_factor, gbuf_2d) if (rootproc) then map_factor_global(1:ix,1:jy) = gbuf_2d(1:ix,1:jy) end if ! temporary increase to dimensions for cori ide=ide+1 jde=jde+1 grid%xp%ide=grid%xp%ide+1 grid%xp%jde=grid%xp%jde+1 call da_patch_to_global(grid, grid%xb%cori, gbuf_2dd) if (rootproc) then cori_global(1:ix+1,1:jy+1) = gbuf_2dd(1:ix+1,1:jy+1) end if ide=ide-1 jde=jde-1 grid%xp%ide=grid%xp%ide-1 grid%xp%jde=grid%xp%jde-1 call da_patch_to_global(grid, grid%xb%landmask, gbuf_2d) if (rootproc) then landmask_global(1:ix,1:jy) = gbuf_2d(1:ix,1:jy) end if #endif if (rootproc) then call da_get_unit(anl_inc_unit) open(unit=anl_inc_unit, file='analysis_increments', form='unformatted') write (unit=anl_inc_unit) ANALYSIS_DATE write (unit=anl_inc_unit) 1, ix, 1, jy, 1, kz ! Map projection information: write (unit=anl_inc_unit) map_projection, coarse_ix, coarse_jy write (unit=anl_inc_unit) & coarse_ds, start_x, start_y, & phic, xlonc, cone_factor, truelat1_3dv, truelat2_3dv, pole, dsm, & psi1, c2, ptop, base_pres, t0, base_lapse, base_temp ! 1d constant fields: write (unit=anl_inc_unit) grid%xb%sigmah, grid%xb%sigmaf #ifdef DM_PARALLEL ! 3d- and 2d-increments: write (unit=anl_inc_unit) u_global, v_global, w_global, p_global, & t_global, q_global, ph_global, mu_global, psfc_global ! 2d-constant fields: write (unit=anl_inc_unit) psac_global, tgrn_global, terr_global, & snow_global, lat_global, lon_global, lanu_global, map_factor_global, & cori_global, landmask_global close(anl_inc_unit) call da_free_unit(anl_inc_unit) #else ! 3d- and 2d-increments: write (unit=anl_inc_unit) grid%xa%u(1:ix+1,1:jy+1,1:kz+1), & grid%xa%v(1:ix+1,1:jy+1,1:kz+1), & grid%xa%w(1:ix+1,1:jy+1,1:kz+1), & grid%xa%p(1:ix+1,1:jy+1,1:kz+1), & grid%xa%t(1:ix+1,1:jy+1,1:kz+1), & q_cgrid(1:ix+1,1:jy+1,1:kz+1), & ph_cgrid(1:ix+1,1:jy+1,1:kz+1), & mu_cgrid(1:ix+1,1:jy+1), & grid%xa%psfc(1:ix+1,1:jy+1) ! .. 2d-constant fields: write (unit=anl_inc_unit) grid%xb%psac(1:ix+1,1:jy+1), & grid%xb%tgrn(1:ix+1,1:jy+1), & grid%xb%terr(1:ix+1,1:jy+1), & grid%xb%snow(1:ix+1,1:jy+1), & grid%xb%lat(1:ix+1,1:jy+1), & grid%xb%lon(1:ix+1,1:jy+1), & grid%xb%lanu(1:ix+1,1:jy+1), & grid%xb%map_factor(1:ix+1,1:jy+1), & grid%xb%cori(1:ix+1,1:jy+1), & grid%xb%landmask(1:ix+1,1:jy+1) close(anl_inc_unit) call da_free_unit(anl_inc_unit) #endif end if if (trace_use) call da_trace_exit("da_write_increments") end subroutine da_write_increments