CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
iccninterp.F90
1
4
9
10 implicit none
11
12 private
13
14 public :: read_cidata, setindxci, ciinterpol
15
16contains
17
19 SUBROUTINE read_cidata (me, master)
20 use machine, only: kind_phys
21 use iccn_def
22 use netcdf
23!--- in/out
24 integer, intent(in) :: me
25 integer, intent(in) :: master
26!--- locals
27 integer :: ncerr
28 integer :: i, n, k, ncid, varid,j,it
29 real(kind=kind_phys), allocatable, dimension(:) :: hyam,hybm
30 real(kind=4), allocatable, dimension(:,:,:) :: ci_ps
31
32 allocate (hyam(kcipl), hybm(kcipl), ci_ps(lonscip,latscip,timeci))
33 allocate (ciplin(lonscip,latscip,kcipl,timeci))
34 allocate (ccnin(lonscip,latscip,kcipl,timeci))
35 allocate (ci_pres(lonscip,latscip,kcipl,timeci))
36 ncerr = nf90_open("cam5_4_143_NAAI_monclimo2.nc", nf90_nowrite, ncid)
37 ncerr = nf90_inq_varid(ncid, "lat", varid)
38 ncerr = nf90_get_var(ncid, varid, ci_lat)
39 ncerr = nf90_inq_varid(ncid, "lon", varid)
40 ncerr = nf90_get_var(ncid, varid, ci_lon)
41 ncerr = nf90_inq_varid(ncid, "PS", varid)
42 ncerr = nf90_get_var(ncid, varid, ci_ps)
43 ncerr = nf90_inq_varid(ncid, "hyam", varid)
44 ncerr = nf90_get_var(ncid, varid, hyam)
45 ncerr = nf90_inq_varid(ncid, "hybm", varid)
46 ncerr = nf90_get_var(ncid, varid, hybm)
47 ncerr = nf90_inq_varid(ncid, "NAAI", varid)
48 ncerr = nf90_get_var(ncid, varid, ciplin)
49 do it = 1,timeci
50 do k=1, kcipl
51 ci_pres(:,:,k,it)=hyam(k)*1.e5+hybm(k)*ci_ps(:,:,it)
52 end do
53 end do
54 ncerr = nf90_close(ncid)
55 ncerr = nf90_open("cam5_4_143_NPCCN_monclimo2.nc", nf90_nowrite, ncid)
56 ncerr = nf90_inq_varid(ncid, "NPCCN", varid)
57 ncerr = nf90_get_var(ncid, varid, ccnin)
58 ncerr = nf90_close(ncid)
59!---
60 deallocate (hyam, hybm, ci_ps)
61 if (me == master) then
62 write(*,*) 'Reading in ICCN data',ci_time
63 endif
64
65 END SUBROUTINE read_cidata
66!
67!**********************************************************************
68!
70 SUBROUTINE setindxci(npts,dlat,jindx1,jindx2,ddy,dlon, &
71 iindx1,iindx2,ddx)
72!
73 USE machine, ONLY: kind_phys
74 USE iccn_def, ONLY: jci => latscip, ci_lat,ici=>lonscip, ci_lon
75!
76 implicit none
77!
78 integer npts, jindx1(npts),jindx2(npts),iindx1(npts),iindx2(npts)
79 real(kind=kind_phys) dlat(npts),ddy(npts),dlon(npts),ddx(npts)
80!
81 integer i,j
82!
83 DO j=1,npts
84 jindx2(j) = jci + 1
85 do i=1,jci
86 if (dlat(j) < ci_lat(i)) then
87 jindx2(j) = i
88 exit
89 endif
90 enddo
91 jindx1(j) = max(jindx2(j)-1,1)
92 jindx2(j) = min(jindx2(j),jci)
93 if (jindx2(j) .ne. jindx1(j)) then
94 ddy(j) = (dlat(j) - ci_lat(jindx1(j))) &
95 / (ci_lat(jindx2(j)) - ci_lat(jindx1(j)))
96 else
97 ddy(j) = 1.0
98 endif
99 !print *,' j=',j,' dlat=',dlat(j),' jindx12=',jindx1(j), &
100 ! jindx2(j),' ci_lat=',ci_lat(jindx1(j)), &
101 ! ci_lat(jindx2(j)),' ddy=',ddy(j)
102 ENDDO
103
104 DO j=1,npts
105 iindx2(j) = ici + 1
106 do i=1,ici
107 if (dlon(j) < ci_lon(i)) then
108 iindx2(j) = i
109 exit
110 endif
111 enddo
112 iindx1(j) = max(iindx2(j)-1,1)
113 iindx2(j) = min(iindx2(j),ici)
114 if (iindx2(j) .ne. iindx1(j)) then
115 ddx(j) = (dlon(j) - ci_lon(iindx1(j))) &
116 / (ci_lon(iindx2(j)) - ci_lon(iindx1(j)))
117 else
118 ddx(j) = 1.0
119 endif
120 !print *,' j=',j,' dlon=',dlon(j),' iindx12=',iindx1(j), &
121 ! iindx2(j),' ci_lon=',ci_lon(iindx1(j)), &
122 ! ci_lon(iindx2(j)),' ddx=',ddx(j)
123 ENDDO
124
125 RETURN
126 END SUBROUTINE setindxci
127!
128!**********************************************************************
129!**********************************************************************
130!
132 SUBROUTINE ciinterpol(me,npts,IDATE,FHOUR,jindx1,jindx2,ddy, &
133 iindx1,iindx2,ddx,lev, prsl, ciplout,ccnout)
134!
135 USE machine, ONLY : kind_phys, kind_dbl_prec
136 use iccn_def
137 implicit none
138 integer i1,i2, iday,j,j1,j2,l,npts,nc,n1,n2,lev,k,i
139 real(kind=kind_phys) fhour,temj, tx1, tx2,temi
140!
141
142 integer jindx1(npts), jindx2(npts),iindx1(npts),iindx2(npts)
143 integer me,idate(4)
144 integer idat(8),jdat(8)
145!
146 real(kind=kind_phys) ddy(npts), ddx(npts),ttt
147 real(kind=kind_phys) ciplout(npts,lev),cipm(npts,kcipl)
148 real(kind=kind_phys) ccnout(npts,lev),ccnpm(npts,kcipl)
149 real(kind=kind_phys) cipres(npts,kcipl), prsl(npts,lev)
150 real(kind=kind_phys) rjday
151 real(kind=kind_dbl_prec) rinc(5)
152 integer jdow, jdoy, jday
153!
154 idat=0
155 idat(1)=idate(4)
156 idat(2)=idate(2)
157 idat(3)=idate(3)
158 idat(5)=idate(1)
159 rinc=0.
160 rinc(2)=fhour
161 CALL w3movdat(rinc,idat,jdat)
162!
163 jdow = 0
164 jdoy = 0
165 jday = 0
166 call w3doxdat(jdat,jdow,jdoy,jday)
167 rjday = jdoy + jdat(5) / 24.
168 IF (rjday .LT. ci_time(1)) rjday = rjday+365.
169!
170 n2 = timeci + 1
171 do j=2,timeci
172 if (rjday .lt. ci_time(j)) then
173 n2 = j
174 exit
175 endif
176 enddo
177 n1 = n2 - 1
178
179!
180!
181 tx1 = (ci_time(n2) - rjday) / (ci_time(n2) - ci_time(n1))
182 if (n2 > timeci) n2 = n2 - timeci
183! if (me .eq. 0) print *,' n1=',n1,' n2=',n2,' rjday=',rjday &
184! ,'ci_time=',ci_time(n1),ci_time(n2), ci_time(timeci+1),tx1
185 tx2 = 1.0 - tx1
186!
187 DO l=1,kcipl
188 DO j=1,npts
189 j1 = jindx1(j)
190 j2 = jindx2(j)
191 temj = 1.0 - ddy(j)
192 i1 = iindx1(j)
193 i2 = iindx2(j)
194 temi = 1.0 - ddx(j)
195 cipm(j,l) = &
196 tx1*(temi*temj*ciplin(i1,j1,l,n1)+ddx(j)*ddy(j)*ciplin(i2,j2,l,n1) &
197 +temi*ddy(j)*ciplin(i1,j2,l,n1)+ddx(j)*temj*ciplin(i2,j1,l,n1)) &
198 + tx2*(temi*temj*ciplin(i1,j1,l,n2)+ddx(j)*ddy(j)*ciplin(i2,j2,l,n2) &
199 +temi*ddy(j)*ciplin(i1,j2,l,n2)+ddx(j)*temj*ciplin(i2,j1,l,n2))
200
201 ccnpm(j,l) = &
202 tx1*(temi*temj*ccnin(i1,j1,l,n1)+ddx(j)*ddy(j)*ccnin(i2,j2,l,n1) &
203 +temi*ddy(j)*ccnin(i1,j2,l,n1)+ddx(j)*temj*ccnin(i2,j1,l,n1)) &
204 + tx2*(temi*temj*ccnin(i1,j1,l,n2)+ddx(j)*ddy(j)*ccnin(i2,j2,l,n2) &
205 +temi*ddy(j)*ccnin(i1,j2,l,n2)+ddx(j)*temj*ccnin(i2,j1,l,n2))
206
207 cipres(j,l) = &
208 tx1*(temi*temj*ci_pres(i1,j1,l,n1)+ddx(j)*ddy(j)*ci_pres(i2,j2,l,n1) &
209 +temi*ddy(j)*ci_pres(i1,j2,l,n1)+ddx(j)*temj*ci_pres(i2,j1,l,n1)) &
210 + tx2*(temi*temj*ci_pres(i1,j1,l,n2)+ddx(j)*ddy(j)*ci_pres(i2,j2,l,n2) &
211 +temi*ddy(j)*ci_pres(i1,j2,l,n2)+ddx(j)*temj*ci_pres(i2,j1,l,n2))
212 ENDDO
213 ENDDO
214
215 DO j=1,npts
216 DO l=1,lev
217 ! noticed input is from top to bottom
218 if(prsl(j,l).ge.cipres(j,kcipl)) then
219 ciplout(j,l)=cipm(j,kcipl)
220 ccnout(j,l)=ccnpm(j,kcipl)
221 else if(prsl(j,l).le.cipres(j,1)) then
222 ciplout(j,l)=cipm(j,1)
223 ccnout(j,l)=ccnpm(j,1)
224 else
225 DO k=kcipl-1,1,-1
226 ! DH* There is no backstop if this condition isn't met,
227 ! i.e. i1 and i2 will have values determined by the
228 ! previous code (line 178) - this leads to crashes in
229 ! debug mode (out of bounds), for example for regression
230 ! test fv3_stretched_nest_debug. For the time being,
231 ! this is 'solved' by simply switching off ICCN
232 ! if MG2/3 are not used (these are the only microphysics
233 ! schemes that use the ICCN data); however, this doesn't
234 ! mean that the code is correct for MG2/3, it just doesn't
235 ! abort if the below condition isn't met, because the code
236 ! is not tested in DEBUG mode. *DH
237 IF(prsl(j,l)>cipres(j,k)) then
238 i1=k
239 i2=min(k+1,kcipl)
240 exit
241 end if
242 end do
243 ciplout(j,l)=cipm(j,i1)+(cipm(j,i2)-cipm(j,i1)) &
244 /(cipres(j,i2)-cipres(j,i1))*(prsl(j,l)-cipres(j,i1))
245 ccnout(j,l)=ccnpm(j,i1)+(ccnpm(j,i2)-ccnpm(j,i1)) &
246 /(cipres(j,i2)-cipres(j,i1))*(prsl(j,l)-cipres(j,i1))
247 end if
248 ENDDO
249 ENDDO
250!
251 RETURN
252 END SUBROUTINE ciinterpol
253
254end module iccninterp
This module defines IN and CCN arrays.
Definition iccn_def.F:6
This module contains subroutines of reading and interplating IN and CCN data.
Definition iccninterp.F90:8