CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
h2ointerp.f90
1
4
9
10 implicit none
11
12 private
13
14 public :: read_h2odata, setindxh2o, h2ointerpol
15
16contains
17
19 subroutine read_h2odata (h2o_phys, me, master)
20 use machine, only: kind_phys
21 use h2o_def
22!--- in/out
23 logical, intent(in) :: h2o_phys
24 integer, intent(in) :: me
25 integer, intent(in) :: master
26!--- locals
27 integer :: i, n, k
28 real(kind=4), allocatable, dimension(:) :: h2o_lat4, h2o_pres4
29 real(kind=4), allocatable, dimension(:) :: h2o_time4, tempin
30
31 if (.not. h2o_phys) then
32 latsh2o = 1
33 levh2o = 1
34 h2o_coeff = 1
35 timeh2o = 1
36
37 return
38 endif
39
40 open(unit=kh2opltc,file='global_h2oprdlos.f77', form='unformatted', convert='big_endian')
41
42!--- read in indices
43!---
44 read (kh2opltc) h2o_coeff, latsh2o, levh2o, timeh2o
45 if (me == master) then
46 write(*,*) 'Reading in h2odata from global_h2oprdlos.f77 '
47 write(*,*) ' h2o_coeff = ', h2o_coeff
48 write(*,*) ' latsh2o = ', latsh2o
49 write(*,*) ' levh2o = ', levh2o
50 write(*,*) ' timeh2o = ', timeh2o
51 endif
52
53!--- read in data
54!--- h2o_lat - latitude of data (-90 to 90)
55!--- h2o_pres - vertical pressure level (mb)
56!--- h2o_time - time coordinate (days)
57!---
58 allocate (h2o_lat(latsh2o), h2o_pres(levh2o),h2o_time(timeh2o+1))
59 allocate (h2o_lat4(latsh2o), h2o_pres4(levh2o),h2o_time4(timeh2o+1))
60 rewind(kh2opltc)
61 read (kh2opltc) h2o_coeff, latsh2o, levh2o, timeh2o, h2o_lat4, h2o_pres4, h2o_time4
62 h2o_pres(:) = h2o_pres4(:)
63!--- convert pressure levels from mb to ln(Pa)
64 h2o_pres(:) = log(100.0*h2o_pres(:))
65 h2o_lat(:) = h2o_lat4(:)
66 h2o_time(:) = h2o_time4(:)
67 deallocate (h2o_lat4, h2o_pres4, h2o_time4)
68
69!--- read in h2oplin which is in order of (lattitudes, water levels, coeff number, time)
70!--- assume latitudes is on a uniform gaussian grid
71!---
72 allocate (tempin(latsh2o))
73 allocate (h2oplin(latsh2o,levh2o,h2o_coeff,timeh2o))
74 DO i=1,timeh2o
75 do n=1,h2o_coeff
76 DO k=1,levh2o
77 READ(kh2opltc) tempin
78 h2oplin(:,k,n,i) = tempin(:)
79 ENDDO
80 enddo
81 ENDDO
82 deallocate (tempin)
83
84 close(kh2opltc)
85
86 end subroutine read_h2odata
87!
88!**********************************************************************
90 subroutine setindxh2o(npts,dlat,jindx1,jindx2,ddy)
91!
92! May 2015 Shrinivas Moorthi - Prepare for H2O interpolation
93!
94 use machine, only: kind_phys
95 use h2o_def, only: jh2o => latsh2o, h2o_lat, h2o_time
96!
97 implicit none
98!
99 integer npts
100 integer, dimension(npts) :: jindx1, jindx2
101 real(kind=kind_phys) :: dlat(npts),ddy(npts)
102!
103 integer i,j,lat
104!
105 do j=1,npts
106 jindx2(j) = jh2o + 1
107 do i=1,jh2o
108 if (dlat(j) < h2o_lat(i)) then
109 jindx2(j) = i
110 exit
111 endif
112 enddo
113 jindx1(j) = max(jindx2(j)-1,1)
114 jindx2(j) = min(jindx2(j),jh2o)
115 if (jindx2(j) /= jindx1(j)) then
116 ddy(j) = (dlat(j) - h2o_lat(jindx1(j))) &
117 / (h2o_lat(jindx2(j)) - h2o_lat(jindx1(j)))
118 else
119 ddy(j) = 1.0
120 endif
121! print *,' j=',j,' dlat=',dlat(j),' jindx12=',jindx1(j), &
122! jindx2(j),' h2o_lat=',h2o_lat(jindx1(j)), &
123! h2o_lat(jindx2(j)),' ddy=',ddy(j)
124 enddo
125
126 return
127 end subroutine setindxh2o
128!
129!**********************************************************************
131 subroutine h2ointerpol(me,npts,idate,fhour,jindx1,jindx2,h2oplout,ddy)
132!
133! May 2015 Shrinivas Moorthi - Prepare for H2O interpolation
134!
135 use machine , only : kind_phys, kind_dbl_prec
136 use h2o_def
137 implicit none
138 integer j,j1,j2,l,npts,nc,n1,n2
139 real(kind=kind_phys) fhour,tem, tx1, tx2
140!
141
142 integer jindx1(npts), jindx2(npts)
143 integer me,idate(4)
144 integer idat(8),jdat(8)
145!
146 real(kind=kind_phys) ddy(npts)
147 real(kind=kind_phys) h2oplout(npts,levh2o,h2o_coeff)
148 real(kind=kind_phys) rjday
149 real(kind=kind_dbl_prec) rinc(5)
150 integer jdow, jdoy, jday
151!
152 idat = 0
153 idat(1) = idate(4)
154 idat(2) = idate(2)
155 idat(3) = idate(3)
156 idat(5) = idate(1)
157 rinc = 0.
158 rinc(2) = fhour
159 CALL w3movdat(rinc,idat,jdat)
160!
161 jdow = 0
162 jdoy = 0
163 jday = 0
164 call w3doxdat(jdat,jdow,jdoy,jday)
165 rjday = jdoy + jdat(5) / 24.
166 if (rjday < h2o_time(1)) rjday = rjday+365.
167!
168 n2 = timeh2o + 1
169 do j=2,timeh2o
170 if (rjday < h2o_time(j)) then
171 n2 = j
172 exit
173 endif
174 enddo
175 n1 = n2 - 1
176!
177! if (me .eq. 0) print *,' n1=',n1,' n2=',n2,' rjday=',rjday
178! &,'h2o_time=',h2o_time(n1),h2o_time(n2)
179!
180 tx1 = (h2o_time(n2) - rjday) / (h2o_time(n2) - h2o_time(n1))
181 tx2 = 1.0 - tx1
182 if (n2 > timeh2o) n2 = n2 - timeh2o
183!
184 do nc=1,h2o_coeff
185 do l=1,levh2o
186 do j=1,npts
187 j1 = jindx1(j)
188 j2 = jindx2(j)
189 tem = 1.0 - ddy(j)
190 h2oplout(j,l,nc) = &
191 tx1*(tem*h2oplin(j1,l,nc,n1)+ddy(j)*h2oplin(j2,l,nc,n1)) &
192 + tx2*(tem*h2oplin(j1,l,nc,n2)+ddy(j)*h2oplin(j2,l,nc,n2))
193 enddo
194 enddo
195 enddo
196!
197 return
198 end subroutine h2ointerpol
199
200end module h2ointerp
This module defines arrays in H2O scheme.
Definition h2o_def.f:6
This module contains subroutines of reading and interpolating h2o coefficients.
Definition h2ointerp.f90:8