var/da/da_setup_structures/da_write_increments.inc

   1 subroutine da_write_increments (grid, q_cgrid, mu_cgrid, ph_cgrid)
   2
   3    !----------------------------------------------------------------------
   4    ! Purpose: Write analysis increments
   5    !----------------------------------------------------------------------
   6
   7    implicit none
   8
   9    type (domain), intent(inout)                         :: grid
  10    real,intent(in) :: q_cgrid(ims:ime,jms:jme,kms:kme)
  11    real,intent(in) :: ph_cgrid(ims:ime,jms:jme,kms:kme)
  12    real,intent(in) :: mu_cgrid(ims:ime,jms:jme)
  13
  14    ! Arrays for write out increments:
  15    integer                                          :: ix, jy, kz
  16 #ifdef DM_PARALLEL
  17    real, dimension(1:grid%xb%mix,1:grid%xb%mjy)               ::     gbuf_2d
  18    real, dimension(1:grid%xb%mix+1,1:grid%xb%mjy+1)           ::     gbuf_2dd
  19    real, dimension(1:grid%xb%mix,1:grid%xb%mjy,1:grid%xb%mkz)      ::     gbuf
  20
  21    real, dimension(1:grid%xb%mix,1:grid%xb%mjy,1:grid%xb%mkz+1)    ::    wgbuf
  22    real, dimension(:,:,:), allocatable :: u_global, v_global, w_global, &
  23                                           p_global, t_global, q_global, &
  24                                          ph_global
  25    real, dimension(:,:)  , allocatable :: mu_global, psfc_global, &
  26                        psac_global, tgrn_global, terr_global, snow_global,&
  27                         lat_global,  lon_global, lanu_global,             &
  28                  map_factor_global, cori_global, landmask_global
  29 #endif
  30
  31    integer :: anl_inc_unit
  32
  33    if (trace_use) call da_trace_entry("da_write_increments")
  34
  35
  36    ! Dimension of the domain:
  37    ix = grid%xb%mix
  38    jy = grid%xb%mjy
  39    kz = grid%xb%mkz
  40
  41 #ifdef DM_PARALLEL
  42
  43    ! 3-d and 2-d increments:
  44
  45    allocate (   p_global (1:ix+1,1:jy+1,1:kz+1))
  46    allocate (   t_global (1:ix+1,1:jy+1,1:kz+1))
  47    allocate (   q_global (1:ix+1,1:jy+1,1:kz+1))
  48    allocate (   u_global (1:ix+1,1:jy+1,1:kz+1))
  49    allocate (   v_global (1:ix+1,1:jy+1,1:kz+1))
  50    allocate (   w_global (1:ix+1,1:jy+1,1:kz+1))
  51    allocate (  ph_global (1:ix+1,1:jy+1,1:kz+1))
  52    allocate (psfc_global (1:ix+1,1:jy+1))
  53    allocate (  mu_global (1:ix+1,1:jy+1))
  54    call da_patch_to_global(grid, grid%xa % p, gbuf)
  55    if (rootproc) then
  56       p_global(1:ix,1:jy,1:kz) = gbuf(1:ix,1:jy,1:kz)
  57    end if
  58    call da_patch_to_global(grid, grid%xa % t, gbuf)
  59    if (rootproc) then
  60       t_global(1:ix,1:jy,1:kz) = gbuf(1:ix,1:jy,1:kz)
  61    end if
  62    call da_patch_to_global(grid, q_cgrid, gbuf)
  63    if (rootproc) then
  64       q_global(1:ix,1:jy,1:kz) = gbuf(1:ix,1:jy,1:kz)
  65    end if
  66    call da_patch_to_global(grid, grid%xa % u, gbuf)
  67    if (rootproc) then
  68       u_global(1:ix,1:jy,1:kz) = gbuf(1:ix,1:jy,1:kz)
  69    end if
  70    call da_patch_to_global(grid, grid%xa % v, gbuf)
  71    if (rootproc) then
  72       v_global(1:ix,1:jy,1:kz) = gbuf(1:ix,1:jy,1:kz)
  73    end if
  74
  75    ! One more level for w and ph:
  76    grid%xp%kde=grid%xp%kde+1
  77    kde=kde+1
  78    call da_patch_to_global(grid, grid%xa % w, wgbuf)
  79    if (rootproc) then
  80       w_global(1:ix,1:jy,1:kz+1) = wgbuf(1:ix,1:jy,1:kz+1)
  81    end if
  82    call da_patch_to_global(grid, ph_cgrid, wgbuf)
  83    if (rootproc) then
  84       ph_global(1:ix,1:jy,1:kz+1) = wgbuf(1:ix,1:jy,1:kz+1)
  85    end if
  86    kde=kde-1
  87    grid%xp%kde=grid%xp%kde-1
  88
  89    call da_patch_to_global(grid, grid%xa % psfc, gbuf_2d)
  90    if (rootproc) then
  91       psfc_global(1:ix,1:jy) = gbuf_2d(1:ix,1:jy)
  92    end if
  93    call da_patch_to_global(grid, mu_cgrid, gbuf_2d)
  94    if (rootproc) then
  95       mu_global(1:ix,1:jy) = gbuf_2d(1:ix,1:jy)
  96    end if
  97
  98    ! 2d constant fields:
  99
 100    allocate (      psac_global (1:ix+1,1:jy+1))
 101    allocate (      tgrn_global (1:ix+1,1:jy+1))
 102    allocate (      terr_global (1:ix+1,1:jy+1))
 103    allocate (      snow_global (1:ix+1,1:jy+1))
 104    allocate (       lat_global (1:ix+1,1:jy+1))
 105    allocate (       lon_global (1:ix+1,1:jy+1))
 106    allocate (      lanu_global (1:ix+1,1:jy+1))
 107    allocate (map_factor_global (1:ix+1,1:jy+1))
 108    allocate (      cori_global (1:ix+1,1:jy+1))
 109    allocate (  landmask_global (1:ix+1,1:jy+1))
 110
 111    call da_patch_to_global(grid, grid%xb%psac, gbuf_2d)
 112    if (rootproc) then
 113       psac_global(1:ix,1:jy) = gbuf_2d(1:ix,1:jy)
 114    end if
 115    call da_patch_to_global(grid, grid%xb%tgrn, gbuf_2d)
 116    if (rootproc) then
 117       tgrn_global(1:ix,1:jy) = gbuf_2d(1:ix,1:jy)
 118    end if
 119    call da_patch_to_global(grid, grid%xb%terr, gbuf_2d)
 120    if (rootproc) then
 121       terr_global(1:ix,1:jy) = gbuf_2d(1:ix,1:jy)
 122    end if
 123    call da_patch_to_global(grid, grid%xb%snow, gbuf_2d)
 124    if (rootproc) then
 125       snow_global(1:ix,1:jy) = gbuf_2d(1:ix,1:jy)
 126    end if
 127    call da_patch_to_global(grid, grid%xb%lat , gbuf_2d)
 128    if (rootproc) then
 129       lat_global(1:ix,1:jy) = gbuf_2d(1:ix,1:jy)
 130    end if
 131    call da_patch_to_global(grid, grid%xb%lon , gbuf_2d)
 132    if (rootproc) then
 133       lon_global(1:ix,1:jy) = gbuf_2d(1:ix,1:jy)
 134    end if
 135    call da_patch_to_global(grid, grid%xb%lanu, gbuf_2d)
 136    if (rootproc) then
 137       lanu_global(1:ix,1:jy) = gbuf_2d(1:ix,1:jy)
 138    end if
 139    call da_patch_to_global(grid, grid%xb%map_factor, gbuf_2d)
 140    if (rootproc) then
 141       map_factor_global(1:ix,1:jy) = gbuf_2d(1:ix,1:jy)
 142    end if
 143
 144    ! temporary increase to dimensions for cori
 145    ide=ide+1
 146    jde=jde+1
 147    grid%xp%ide=grid%xp%ide+1
 148    grid%xp%jde=grid%xp%jde+1
 149    call da_patch_to_global(grid, grid%xb%cori, gbuf_2dd)
 150    if (rootproc) then
 151       cori_global(1:ix+1,1:jy+1) = gbuf_2dd(1:ix+1,1:jy+1)
 152    end if
 153    ide=ide-1
 154    jde=jde-1
 155    grid%xp%ide=grid%xp%ide-1
 156    grid%xp%jde=grid%xp%jde-1
 157
 158    call da_patch_to_global(grid, grid%xb%landmask, gbuf_2d)
 159    if (rootproc) then
 160       landmask_global(1:ix,1:jy) = gbuf_2d(1:ix,1:jy)
 161    end if
 162
 163 #endif
 164
 165    if (rootproc) then
 166       call da_get_unit(anl_inc_unit)
 167       open(unit=anl_inc_unit, file='analysis_increments', form='unformatted')
 168
 169       write (unit=anl_inc_unit) ANALYSIS_DATE
 170
 171       write (unit=anl_inc_unit) 1, ix, 1, jy, 1, kz
 172
 173       ! Map projection information:
 174       write (unit=anl_inc_unit) map_projection, coarse_ix, coarse_jy
 175       write (unit=anl_inc_unit) &
 176          coarse_ds, start_x, start_y, &
 177          phic, xlonc, cone_factor, truelat1_3dv, truelat2_3dv, pole, dsm,   &
 178          psi1, c2, ptop, base_pres, t0, base_lapse, base_temp
 179
 180       ! 1d constant fields:
 181
 182       write (unit=anl_inc_unit) grid%xb%sigmah, grid%xb%sigmaf
 183
 184 #ifdef DM_PARALLEL
 185
 186       ! 3d- and 2d-increments:
 187       write (unit=anl_inc_unit) u_global, v_global, w_global, p_global, &
 188          t_global, q_global, ph_global, mu_global, psfc_global
 189
 190       ! 2d-constant fields:
 191       write (unit=anl_inc_unit) psac_global, tgrn_global, terr_global, &
 192          snow_global, lat_global, lon_global, lanu_global, map_factor_global, &
 193          cori_global, landmask_global
 194       close(anl_inc_unit)
 195       call da_free_unit(anl_inc_unit)
 196 #else
 197
 198       ! 3d- and 2d-increments:
 199       write (unit=anl_inc_unit) grid%xa%u(1:ix+1,1:jy+1,1:kz+1), &
 200                     grid%xa%v(1:ix+1,1:jy+1,1:kz+1), &
 201                     grid%xa%w(1:ix+1,1:jy+1,1:kz+1), &
 202                     grid%xa%p(1:ix+1,1:jy+1,1:kz+1), &
 203                     grid%xa%t(1:ix+1,1:jy+1,1:kz+1), &
 204                     q_cgrid(1:ix+1,1:jy+1,1:kz+1), &
 205                     ph_cgrid(1:ix+1,1:jy+1,1:kz+1), &
 206                     mu_cgrid(1:ix+1,1:jy+1), &
 207                     grid%xa%psfc(1:ix+1,1:jy+1)
 208
 209       !    .. 2d-constant fields:
 210       write (unit=anl_inc_unit) grid%xb%psac(1:ix+1,1:jy+1), &
 211                     grid%xb%tgrn(1:ix+1,1:jy+1), &
 212                     grid%xb%terr(1:ix+1,1:jy+1), &
 213                     grid%xb%snow(1:ix+1,1:jy+1), &
 214                     grid%xb%lat(1:ix+1,1:jy+1), &
 215                     grid%xb%lon(1:ix+1,1:jy+1), &
 216                     grid%xb%lanu(1:ix+1,1:jy+1), &
 217                     grid%xb%map_factor(1:ix+1,1:jy+1), &
 218                     grid%xb%cori(1:ix+1,1:jy+1), &
 219                     grid%xb%landmask(1:ix+1,1:jy+1)
 220       close(anl_inc_unit)
 221       call da_free_unit(anl_inc_unit)
 222 #endif
 223
 224    end if
 225
 226    if (trace_use) call da_trace_exit("da_write_increments")
 227
 228 end subroutine da_write_increments
 229
 230