da_write_increments.inc
References to this file elsewhere.
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