da_write_kma_increments.inc
References to this file elsewhere.
1 subroutine da_write_kma_increments(xp, xb, xa)
2
3 !---------------------------------------------------------------------------
4 ! Purpose: Gathers KMA analysis increments and writes
5 ! on i"anl_inc_unit" unit
6 !---------------------------------------------------------------------------
7
8 implicit none
9
10 type(xb_type),intent(in) :: xb
11 type(x_type),intent(in) :: xa
12 type(xpose_type),intent(in) :: xp
13
14 ! Arrays for write out increments:
15 integer :: ix, jy, kz
16
17 #ifdef DM_PARALLEL
18 real, dimension(1:xb%mix,1:xb%mjy) :: gbuf_2d
19 real, dimension(1:xb%mix,1:xb%mjy,1:xb%mkz) :: gbuf
20 real, dimension(:,:) , allocatable :: psfc_g
21 real, dimension(:,:,:), allocatable :: u_g, v_g, t_g, q_g, p_g
22 #endif
23
24 integer :: i, j, k,anl_inc_unit
25
26 if (trace_use) call da_trace_entry("da_write_kma_increments")
27
28 ! Dimension of the domain:
29 ix = xb%mix
30 jy = xb%mjy
31 kz = xb%mkz
32
33 #ifdef DM_PARALLEL
34
35 ! 3-d and 2-d increments:
36
37 allocate (psfc_g (1:ix,1:jy))
38 allocate ( u_g (1:ix,1:jy,1:kz))
39 allocate ( v_g (1:ix,1:jy,1:kz))
40 allocate ( t_g (1:ix,1:jy,1:kz))
41 allocate ( q_g (1:ix,1:jy,1:kz))
42 allocate ( p_g (1:ix,1:jy,1:kz))
43
44 call da_local_to_global(xp, xa % psfc, gbuf_2d, 2)
45 if (rootproc) then
46 psfc_g(1:ix,1:jy) = gbuf_2d(1:ix,1:jy)
47 end if
48
49 call da_local_to_global(xp, xa % u, gbuf, 3)
50 if (rootproc) then
51 u_g(1:ix,1:jy,1:kz) = gbuf(1:ix,1:jy,1:kz)
52 end if
53
54 call da_local_to_global(xp, xa % v, gbuf, 3)
55 if (rootproc) then
56 v_g(1:ix,1:jy,1:kz) = gbuf(1:ix,1:jy,1:kz)
57 end if
58
59 call da_local_to_global(xp, xa % t, gbuf, 3)
60 if (rootproc) then
61 t_g(1:ix,1:jy,1:kz) = gbuf(1:ix,1:jy,1:kz)
62 end if
63
64 call da_local_to_global(xp, xa % q, gbuf, 3)
65 if (rootproc) then
66 q_g(1:ix,1:jy,1:kz) = gbuf(1:ix,1:jy,1:kz)
67 end if
68
69 call da_local_to_global(xp, xa % p, gbuf, 3)
70 if (rootproc) then
71 p_g(1:ix,1:jy,1:kz) = gbuf(1:ix,1:jy,1:kz)
72 end if
73 #endif
74
75 if (rootproc) then
76 ! 3d- and 2d-increments:
77
78 call da_get_unit(anl_inc_unit)
79 open(unit=anl_inc_unit,file="analysis_increments_kma",status="replace", &
80 form="unformatted")
81 #ifdef DM_PARALLEL
82 write(anl_inc_unit) &
83 ((psfc_g(i,j),i=xb%ids,xb%ide),j=xb%jds,xb%jde)
84 write(anl_inc_unit) &
85 (((u_g(i,j,k),i=xb%ids,xb%ide),j=xb%ids,xb%jde),k=xb%kds,xb%kde)
86 write(anl_inc_unit) &
87 (((v_g(i,j,k),i=xb%ids,xb%ide),j=xb%ids,xb%jde),k=xb%kds,xb%kde)
88 write(anl_inc_unit) &
89 (((t_g(i,j,k),i=xb%ids,xb%ide),j=xb%ids,xb%jde),k=xb%kds,xb%kde)
90 write(anl_inc_unit) &
91 (((q_g(i,j,k),i=xb%ids,xb%ide),j=xb%ids,xb%jde),k=xb%kds,xb%kde)
92 write(anl_inc_unit) &
93 (((p_g(i,j,k),i=xb%ids,xb%ide),j=xb%ids,xb%jde),k=xb%kds,xb%kde)
94
95 write(unit=stdout,fmt='(10e15.7)') &
96 ((psfc_g(i,j),i=xb%ids,xb%ide),j=xb%jds,xb%jde)
97 write(unit=stdout,fmt='(10e15.7)') &
98 (((u_g(i,j,k),i=xb%ids,xb%ide),j=xb%ids,xb%jde),k=xb%kds,xb%kde)
99 write(unit=stdout,fmt='(10e15.7)') &
100 (((v_g(i,j,k),i=xb%ids,xb%ide),j=xb%ids,xb%jde),k=xb%kds,xb%kde)
101 write(unit=stdout,fmt='(10e15.7)') &
102 (((t_g(i,j,k),i=xb%ids,xb%ide),j=xb%ids,xb%jde),k=xb%kds,xb%kde)
103 write(unit=stdout,fmt='(10e15.7)') &
104 (((q_g(i,j,k),i=xb%ids,xb%ide),j=xb%ids,xb%jde),k=xb%kds,xb%kde)
105 write(unit=stdout,fmt='(10e15.7)') &
106 (((p_g(i,j,k),i=xb%ids,xb%ide),j=xb%ids,xb%jde),k=xb%kds,xb%kde)
107 #else
108 write(anl_inc_unit) &
109 ((xa%psfc(i,j),i=xb%ids,xb%ide),j=xb%jds,xb%jde)
110 write(anl_inc_unit) &
111 (((xa%u(i,j,k),i=xb%ids,xb%ide),j=xb%jds,xb%jde),k=xb%kds,xb%kde)
112 write(anl_inc_unit) &
113 (((xa%v(i,j,k),i=xb%ids,xb%ide),j=xb%jds,xb%jde),k=xb%kds,xb%kde)
114 write(anl_inc_unit) &
115 (((xa%t(i,j,k),i=xb%ids,xb%ide),j=xb%jds,xb%jde),k=xb%kds,xb%kde)
116 write(anl_inc_unit) &
117 (((xa%q(i,j,k),i=xb%ids,xb%ide),j=xb%jds,xb%jde),k=xb%kds,xb%kde)
118 write(anl_inc_unit) &
119 (((xa%p(i,j,k),i=xb%ids,xb%ide),j=xb%jds,xb%jde),k=xb%kds,xb%kde)
120
121 write(unit=stdout,fmt='(10e15.7)') &
122 ((xa%psfc(i,j),i=xb%ids,xb%ide),j=xb%jds,xb%jde)
123 write(unit=stdout,fmt='(10e15.7)') &
124 (((xa%u(i,j,k),i=xb%ids,xb%ide),j=xb%jds,xb%jde),k=xb%kds,xb%kde)
125 write(unit=stdout,fmt='(10e15.7)') &
126 (((xa%v(i,j,k),i=xb%ids,xb%ide),j=xb%jds,xb%jde),k=xb%kds,xb%kde)
127 write(unit=stdout,fmt='(10e15.7)') &
128 (((xa%t(i,j,k),i=xb%ids,xb%ide),j=xb%jds,xb%jde),k=xb%kds,xb%kde)
129 write(unit=stdout,fmt='(10e15.7)') &
130 (((xa%q(i,j,k),i=xb%ids,xb%ide),j=xb%jds,xb%jde),k=xb%kds,xb%kde)
131 write(unit=stdout,fmt='(10e15.7)') &
132 (((xa%p(i,j,k),i=xb%ids,xb%ide),j=xb%jds,xb%jde),k=xb%kds,xb%kde)
133 #endif
134 close(anl_inc_unit)
135 call da_free_unit(anl_inc_unit)
136
137 end if
138
139 if (trace_use) call da_trace_exit("da_write_kma_increments")
140
141 end subroutine da_write_kma_increments
142
143