19 SUBROUTINE read_cidata (me, master)
24 integer,
intent(in) :: me
25 integer,
intent(in) :: master
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
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)
51 ci_pres(:,:,k,it)=hyam(k)*1.e5+hybm(k)*ci_ps(:,:,it)
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)
60 deallocate (hyam, hybm, ci_ps)
61 if (me == master)
then
62 write(*,*)
'Reading in ICCN data',ci_time
70 SUBROUTINE setindxci(npts,dlat,jindx1,jindx2,ddy,dlon, &
74 USE iccn_def,
ONLY: jci => latscip, ci_lat,ici=>lonscip, ci_lon
78 integer npts, jindx1(npts),jindx2(npts),iindx1(npts),iindx2(npts)
79 real(kind=kind_phys) dlat(npts),ddy(npts),dlon(npts),ddx(npts)
86 if (dlat(j) < ci_lat(i))
then
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)))
107 if (dlon(j) < ci_lon(i))
then
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)))
132 SUBROUTINE ciinterpol(me,npts,IDATE,FHOUR,jindx1,jindx2,ddy, &
133 iindx1,iindx2,ddx,lev, prsl, ciplout,ccnout)
135 USE machine,
ONLY : kind_phys, kind_dbl_prec
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
142 integer jindx1(npts), jindx2(npts),iindx1(npts),iindx2(npts)
144 integer idat(8),jdat(8)
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
161 CALL w3movdat(rinc,idat,jdat)
166 call w3doxdat(jdat,jdow,jdoy,jday)
167 rjday = jdoy + jdat(5) / 24.
168 IF (rjday .LT. ci_time(1)) rjday = rjday+365.
172 if (rjday .lt. ci_time(j))
then
181 tx1 = (ci_time(n2) - rjday) / (ci_time(n2) - ci_time(n1))
182 if (n2 > timeci) n2 = n2 - timeci
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))
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))
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))
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)
237 IF(prsl(j,l)>cipres(j,k))
then
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))