File: C:\NOAA\NEMS_11731\src\atmos\gfs\dyn\grid_to_spect_inp_1.f

1           subroutine grid_to_spect_inp_1
2          &    (psg,uug,vvg,ttg,rqg,
3          &     trie_ps,trio_ps,
4          &     trie_di,trio_di,trie_ze,trio_ze,
5          &     trie_te,trio_te,trie_rq,trio_rq,
6          &     ls_node,ls_nodes,max_ls_nodes,
7          &     lats_nodes_a,global_lats_a,lonsperlat,
8          &     epse,epso, plnew_a,plnow_a,pwat,ptot,ptrc)
9     !!
10     !! hmhj - this routine do spectral to grid transform 
11     !!        from gfsio read in field, to model fields
12     !! input psg,uug,vvg,ttg,rqg (mapping wind, temp)
13     !! output psg,uug,vvg,ttg,rqg in model values (mapping wind, enthalpy)
14     !! aug 2010      sarah lu, modified to compute tracer global sum
15     !!
16           use gfs_dyn_resol_def
17           use gfs_dyn_layout1
18           use gfs_dyn_gg_def
19           use gfs_dyn_vert_def
20           use gfs_dyn_date_def
21           use namelist_dynamics_def
22           use gfs_dyn_coordinate_def 
23           use gfs_dyn_tracer_const
24           use gfs_dyn_tracer_config, only: glbsum                     !glbsum
25           use gfs_dyn_physcons, fv => con_fvirt, rerth => con_rerth,
26          &              grav => con_g,  cp => con_cp , rd => con_rd
27           implicit none
28     !!
29           real(kind=kind_evod) psg(lonf,lats_node_a)
30           real(kind=kind_evod) uug(lonf,lats_node_a,levs)
31           real(kind=kind_evod) vvg(lonf,lats_node_a,levs)
32           real(kind=kind_evod) ttg(lonf,lats_node_a,levs)
33           real(kind=kind_evod) rqg(lonf,lats_node_a,levh)
34     !
35           REAL(KIND=KIND_GRID) pwat   (lonf,lats_node_a)
36           REAL(KIND=KIND_GRID) ptot   (lonf,lats_node_a)
37           REAL(KIND=KIND_GRID) ptrc   (lonf,lats_node_a,ntrac)        !glbsum
38           REAL(KIND=KIND_GRID) work   (lonf)
39           REAL(KIND=KIND_GRID) tki    (lonf,levp1)
40           REAL(KIND=KIND_GRID) prsi   (lonf,levp1)
41     
42           real(kind=kind_evod)  tkrt0
43           real(kind=kind_evod), parameter :: rkappa = cp / rd
44     !
45           real(kind=kind_evod) trie_ps(len_trie_ls,2)
46           real(kind=kind_evod) trio_ps(len_trio_ls,2)
47           real(kind=kind_evod) trie_di(len_trie_ls,2,levs)
48           real(kind=kind_evod) trio_di(len_trio_ls,2,levs)
49           real(kind=kind_evod) trie_ze(len_trie_ls,2,levs)
50           real(kind=kind_evod) trio_ze(len_trio_ls,2,levs)
51           real(kind=kind_evod) trie_te(len_trie_ls,2,levs)
52           real(kind=kind_evod) trio_te(len_trio_ls,2,levs)
53           real(kind=kind_evod) trie_rq(len_trie_ls,2,levh)
54           real(kind=kind_evod) trio_rq(len_trio_ls,2,levh)
55     !
56     !!!!  integer, parameter :: lota = 3*levs+1*levh+1 
57     !
58           real(kind=kind_evod) trie_ls(len_trie_ls,2,lota)
59           real(kind=kind_evod) trio_ls(len_trio_ls,2,lota)
60     !!
61           real(kind=kind_evod) for_gr_a_1(lonfx*(lota),lats_dim_a)
62           real(kind=kind_evod) for_gr_a_2(lonfx*(lota),lats_dim_a)
63     !
64           integer              ls_node(ls_dim,3)
65           integer              ls_nodes(ls_dim,nodes)
66           integer              max_ls_nodes(nodes)
67           integer              lats_nodes_a(nodes)
68           integer              global_lats_a(latg)
69           integer                 lonsperlat(latg)
70           integer dimg
71     !
72           real(kind=kind_evod)  epse(len_trie_ls)
73           real(kind=kind_evod)  epso(len_trio_ls)
74     !
75           real(kind=kind_evod)   plnew_a(len_trie_ls,latg2)
76           real(kind=kind_evod)   plnow_a(len_trio_ls,latg2)
77     !
78           real(kind=kind_evod)   tfac(lonf,levs), sumq(lonf,levs), rcs2
79     !
80           integer              i,j,k,kk, nn, nnl
81           integer              l,lan,lat
82           integer              lon_dim,lons_lat
83     !
84           integer              locl,n
85           integer              indev
86           integer              indod
87           integer              indev1,indev2
88           integer              indod1,indod2
89           INTEGER              INDLSEV,JBASEV
90           INTEGER              INDLSOD,JBASOD
91     !
92           logical 	lslag
93           logical , parameter :: repro = .false.
94     !
95     
96           real(kind=kind_evod), parameter :: one=1.0D0, pa2cb=0.001D0
97     !
98     !timers______________________________________________________---
99           real*8 rtc ,timer1,timer2
100     !timers______________________________________________________---
101     !
102     !
103           real(kind=kind_evod), parameter :: cons_0=0.0D0, cons_24=24.0D0
104          &,                                cons_99=99.0D0, cons_1p0d9=1.0D9
105          &,                                  qmin=1.0D-10
106     !
107           real(kind=kind_evod) ga2, tem
108     !
109           INCLUDE 'function2'
110     
111     !
112     !--------------------------------------------------------------------
113     !
114           lslag   = .false.
115     !
116           trie_ls = 0.0D0
117           trio_ls = 0.0D0
118     !
119     !--------------------------------------------------------------------
120           do lan=1,lats_node_a
121             lon_dim = lon_dims_a(lan)
122             lat = global_lats_a(ipt_lats_node_a-1+lan)
123             lons_lat = lonsperlat(lat)
124             rcs2     = rcs2_a(min(lat,latg-lat+1))
125     !
126             if (thermodyn_id == 3) then
127               do k=1,levs
128                 do i=1,lons_lat
129                   tfac(i,k) = 0.0D0
130                   sumq(i,k) = 0.0D0
131                 enddo
132               enddo
133               do nn=1,ntrac
134                 nnl = (nn-1)*levs
135                 if (cpi(nn) .ne. 0.0) then
136                   do k=1,levs
137                     do i=1,lons_lat
138                       sumq(i,k) = sumq(i,k) + rqg(i,lan,nnl+k)
139                       tfac(i,k) = tfac(i,k) + cpi(nn)*rqg(i,lan,nnl+k)
140                     enddo
141                   enddo
142                 endif
143               enddo
144               do k=1,levs
145                 do i=1,lons_lat
146                   tfac(i,k) = (one-sumq(i,k))*cpi(0) + tfac(i,k)
147                 enddo
148               enddo
149             else
150               do k=1,levs
151                 do i=1,lons_lat
152                   tfac(i,k) = one + fv*max(rqg(i,lan,k),qmin) 
153                 enddo
154               enddo
155             endif
156     
157             do k=1,levs
158               do i=1,lons_lat
159                 uug(i,lan,k) = uug(i,lan,k) * coslat_a(lat)
160                 vvg(i,lan,k) = vvg(i,lan,k) * coslat_a(lat)
161                 ttg(i,lan,k) = ttg(i,lan,k) * tfac(i,k)
162                 for_gr_a_2(i+(kat+k-2)*lon_dim,lan) = ttg(i,lan,k)
163                 for_gr_a_2(i+(kau+k-2)*lon_dim,lan) = uug(i,lan,k) * rcs2
164                 for_gr_a_2(i+(kav+k-2)*lon_dim,lan) = vvg(i,lan,k) * rcs2
165               enddo
166             enddo
167             do k=1,levh
168               do i=1,lons_lat
169                 for_gr_a_2(i+(kar+k-2)*lon_dim,lan)=rqg(i,lan,k)
170               enddo
171             enddo
172             do i=1,lons_lat
173               ptot(i,lan) = psg(i,lan) * pa2cb
174             enddo
175             if (gen_coord_hybrid) then   ! Ps is the prognostic variable
176               do i=1,lons_lat
177                 psg(i,lan) = psg(i,lan) * pa2cb
178               enddo
179             else                         ! ln(Ps) is the prognostic variable
180               do i=1,lons_lat
181                 psg(i,lan) = log(psg(i,lan)*pa2cb)
182               enddo
183             endif
184             do i=1,lons_lat
185               for_gr_a_2(i+(kaps-1)*lon_dim,lan) = psg(i,lan)
186             enddo
187     !
188     ! get pressure at interfaces for pwat 
189             if (gen_coord_hybrid) then  
190               tki = 0.0D0
191               do k=2,levs
192                 do i=1,lons_lat
193                   tkrt0 = (ttg(i,lan,k-1)+ttg(i,lan,k))
194          &                           /(thref(k-1)+thref(k))
195                   tki (i,k)=ck5(k)*tkrt0**rkappa
196                 enddo
197               enddo
198               do k=1,levp1
199                 do i=1,lons_lat
200                   prsi(i,k)  = ak5(k)+bk5(k)*psg(i,lan)+tki(i,k) 
201                 enddo
202               enddo
203             else if (hybrid) then
204               do k=1,levp1
205                 kk=levp1+1-k
206                 do i=1,lons_lat
207                   prsi(i,k)  = ak5(kk)+bk5(kk)*ptot(i,lan)
208                 enddo
209               enddo
210             else
211               do k=1,levp1
212                 do i=1,lons_lat
213                   prsi(i,k)  = si(k)*ptot(i,lan)
214                 enddo
215               enddo
216             endif                      
217     !
218     ! get pwat (total vertical integrated water)
219             do i=1,lons_lat
220               pwat(i,lan) = 0.0D0
221             enddo
222             do k=1,levs
223               do i=1,lons_lat
224                 work(i) = 0.0D0
225               enddo
226               if( ncld.gt.0 ) then
227                 do nn=ntcw,ntcw+ncld-1
228                   nnl = (nn-1)*levs
229                   do i=1,lons_lat
230                     work(i) = work(i) + rqg(i,lan,nnl+k)
231                   enddo
232                 enddo
233               endif
234               do i=1,lons_lat
235                 pwat(i,lan) = pwat(i,lan) + (prsi(i,k)-prsi(i,k+1))
236          &                                * (rqg(i,lan,k) + work(i))
237               enddo
238             enddo
239     !
240     ! get ptrc (tracer global sum)                                   !glbsum
241     !
242             if( glbsum ) then                                        !glbsum
243               do nn = 1, ntrac                                       !glbsum
244                 nnl = (nn-1)*levs                                    !glbsum
245                 do i=1,lons_lat                                      !glbsum
246                  ptrc(i,lan,nn) = 0.0D0                              !glbsum
247                  do k=1,levs                                         !glbsum
248                    ptrc(i,lan,nn) = ptrc(i,lan,nn) +                 !glbsum
249          &         (prsi(i,k)-prsi(i,k+1))*rqg(i,lan,nnl+k)          !glbsum
250                  enddo                                               !glbsum
251                 enddo                                                !glbsum
252               enddo                                                  !glbsum
253             endif                                                    !glbsum
254     
255     !
256           enddo
257     !
258     ! =======================================================================
259           do lan=1,lats_node_a
260     !
261              lon_dim = lon_dims_a(lan)
262     !
263              lat = global_lats_a(ipt_lats_node_a-1+lan)
264              lons_lat = lonsperlat(lat)
265     
266              call grid2four_thread(for_gr_a_2(1,lan),for_gr_a_1(1,lan),
267          &                  lon_dim,lons_lat,lonfx,lota)
268     !
269           enddo
270     !
271           dimg=0
272           call four2fln(lslag,lats_dim_a,lota,lota,for_gr_a_1,
273          x              ls_nodes,max_ls_nodes,
274          x              lats_nodes_a,global_lats_a,lon_dims_a,
275          x              lats_node_a,ipt_lats_node_a,dimg,
276          x              lat1s_a,lonfx,latg,latg2,
277          x              trie_ls(1,1,1), trio_ls(1,1,1),
278          x              plnew_a, plnow_a,
279          x              ls_node)
280     !
281     !
282           trie_di = 0.0D0
283           trio_di = 0.0D0
284           trie_ze = 0.0D0
285           trio_ze = 0.0D0
286     !
287     !$omp parallel do shared(trie_ls,trio_ls)
288     !$omp+shared(trie_di,trio_di,trie_ze,trio_ze,trie_te,trio_te)
289     !$omp+shared(kau,kav,kat,epse,epso,ls_node)
290     !$omp+private(k)
291           do k=1,levs
292              call uveodz(trie_ls(1,1,kau+k-1), trio_ls(1,1,kav+k-1),
293          x               trie_di(1,1,k),       trio_ze(1,1,k),
294          x               epse,epso,ls_node)
295     !
296              call uvoedz(trio_ls(1,1,kau+k-1), trie_ls(1,1,kav+k-1),
297          x               trio_di(1,1,k),       trie_ze(1,1,k),
298          x               epse,epso,ls_node)
299             trie_te(:,:,k)=trie_ls(:,:,kat+k-1)
300             trio_te(:,:,k)=trio_ls(:,:,kat+k-1)
301           enddo
302           do k=1,levh
303             trie_rq(:,:,k)=trie_ls(:,:,kar+k-1)
304             trio_rq(:,:,k)=trio_ls(:,:,kar+k-1)
305           enddo
306           trie_ps(:,:)=trie_ls(:,:,kaps)
307           trio_ps(:,:)=trio_ls(:,:,kaps)
308     
309           END SUBROUTINE grid_to_spect_inp_1
310