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