program uv_flux c c*********************************************************************** c c author c colin seftor c c date c 28-mar-2001 c c purpose c calculate surface uv flux for 305, 310, 325, and 380 nm c and calculate erythemal flux c c method c performs table lookups on pre-computed flux values using c the TOMS Version 7 ozone and temperature profiles and c interpolation techniques (see TOMS V7 users guide for c more detail) c c variables c name type description c ---- ---- ----------- c cofs r(4,4) lagrange coefficients c cthet0 r(4) used in lagrange coef calc c ctheta r(4) used in lagrange coef calc c daybnd r(12) middle of the month day number c densol r(4,7) used in lagrange coef calc c denscn r(4,3) used in lagrange coef calc c lat real latitude c long real longitude c latind int(2) latitude index (used in profile mixing fraction) c minref real minimum surface reflectivity c p1 real phase function used in reflectivity int c pr real phase function used in reflectivity int c snow int snow/ice flag c surf int surface type c uvf r(4) uv surface flux for 305, 310, 325, 380 nm c eryflux r(4) erythemal flux c c include files c minref.inc minimum reflectivity data base c const.inc constant common block c tab.inc 380 nm radiance table c uvflux.inc uv flux radiance table c c*********************************************************************** c include 'minref.inc' include 'const.inc' include 'tab.inc' include 'uvflux.inc' real*8 lat, long, uvf(4), cluvf(4), eryflux, eryclear real*8 daybnd(12), p1, pr, minref real*8 cofs(16), cthet0(4), ctheta(4) real*8 densol(4,7), denscn(4,3) integer surf, algflg, snow, latind(2), errflg logical success data daybnd/15.,46.,74.,105.,135.,166.,196.,227.,258.,288., 1 319.,349./ c -- wav(1) = 305.088716 wav(2) = 310.089969 wav(3) = 324.093502 wav(4) = 380.107892 open (20,file='ery.dat',status='old') do i = 1,22 read (20,*) ww,ery(i) wav(i+4) = ww enddo close (20) c -- open and read in values open (21,file='minref.ed',status='old',form='unformatted') read (21) mref close (21) open (22,file='TABLE.NEWFLUX',status='old',form = 'unformatted') read (22) uvflux close (22) open (23,file='TABLE.N7',status='old',form = 'unformatted') read (23) logi0n, z1i0n, z2i0n, ti0n, sbn read (23) li0rn, z1i0rn, z2i0rn, ti0rn, sbrn close (23) c -- initialize ccf routines call init(380,success) c -- angle conversions rtod=180./3.14159 dtor=3.14159/180. c -- log pressures p04log = dlog(0.4) p10log = dlog(1.0) c -- calculate values needed in table interpolation do j = 1,7 do k = j,j+3 xdenom = 1.0 do i = j,j+3 if (i.ne.k) xdenom = xdenom * (xzlog(k) - xzlog(i)) enddo densol(k-j+1,j) = xdenom enddo enddo do j = 1, 3 do k = j, j + 3 xdenom = 1.0 do i = j,j+3 if (i.ne.k) xdenom = xdenom * (xlog(k) - xlog(i)) enddo denscn(k-j+1,j) = xdenom enddo enddo c -- open test file open (24,file='global.input',status='old') c open (24,file='test.input',status='old') c -- start reading values from test file c do i = 1,35 c read (24,*) c enddo read (24,*) write (6,*) '' do irec = 1,300000 read (24,*,end=900) iyr,iday,itime,isamp,satza,lat,long, 1 pteran,pcloud,sza,phi,ozone,v7ref,res,soi,xnval,algflg, 2 errflg,surf,snow,psi c -- satellite zenith angle c satza = rtod*asin(1.15*sin(dtor*3.*abs(isamp-18))) c -- convert 380 n value to radiance rad = 10.**(-.01*xnval) c -- day of year doy = float(iday) c -- sun-earth distance phase = 0.0172024238D0 * (dfloat(doy)-3.4) r = 1.0 - 0.0167*dcos(phase) r2 = r*r c -- minimum uv reflectivity if (doy.le.daybnd(1)) then frac = doy/31. bmonth = 12 emonth = 1 else if (doy.ge.daybnd(12)) then frac = (doy-daybnd(12))/31. bmonth = 12 emonth = 1 else do j = 1,11 if (doy.ge.daybnd(j) .and. doy.lt.daybnd(j+1)) then frac = (doy - daybnd(j)) / (daybnd(j+1)-daybnd(j)) bmonth = j emonth = j+1 goto 10 endif enddo 10 continue endif mlat = lat + 91. if (mlat.gt.180) mlat = 180 mlon = (long + 181.)/1.25 if (mlon.gt.288) mlon = 288 if (mref(mlon,mlat,bmonth).ne.990. .and. 1 mref(mlon,mlat,emonth).ne.990.) then minref = ((1.-frac)*mref(mlon,mlat,bmonth) + 1 frac*mref(mlon,mlat,emonth)) / 100. else if (mref(mlon,mlat,bmonth).eq.990. .and. 1 mref(mlon,mlat,emonth).ne.990.) then minref = mref(mlon,mlat,emonth) / 100. else if (mref(mlon,mlat,bmonth).ne.990. .and. 1 mref(mlon,mlat,emonth).eq.990.) then minref = mref(mlon,mlat,bmonth) / 100. else if (abs(lat) .gt. 65) then minref = 0.70 else minref = 0.03 endif c -- table position do i=1,10 indsol=i if (xzlog(i).ge.dlog(1./cosd(sza))) go to 20 enddo 20 continue indsol=indsol-2 if(indsol.lt.1) indsol=1 if(indsol.gt.7) indsol=7 do i=1,6 indscn=i if (xlog(i).ge.dlog(1./cosd(satza))) go to 30 enddo 30 continue indscn=indscn-2 if(indscn.lt.1) indscn=1 if(indscn.gt.3) indscn=3 iofset=(indsol-1)*6+(indscn-1) c -- lagrange coefficients indmax=indsol+3 j=1 xzlog1 = 0. xlog1 = 0. if (cosd(sza).ne.0.) xzlog1 = dlog(1.0/dcosd(sza)) if (cosd(satza).ne.0) xlog1 = dlog(1.0/dcosd(satza)) do k=indsol,indmax xnom=1.0 do i=indsol,indmax if (i.ne.k) xnom = (xzlog1-xzlog(i))*xnom enddo cthet0(j)=xnom/densol(j,indsol) j=j+1 enddo j=1 indmax=indscn+3 do k=indscn,indmax xnom=1.0 do i=indscn,indmax if (i.ne.k) xnom=(xlog1-xlog(i))*xnom enddo ctheta(j)=xnom/denscn(j,indscn) j=j+1 enddo l = 1 do i=1,4 do k=1,4 cofs(l)=ctheta(k)*cthet0(i) l=l+1 enddo enddo c -- phase factors p1 = -.375*cosd(sza)*sind(sza)*sind(satza) pr = 2.0*p1/(3.0*cosd(satza)*cosd(sza)*cosd(sza)) c -- aerosol screen c psi = acos(cosd(sza)*cosd(satza) + c 1 sind(sza)*sind(satza)*cosd(phi)) c psi = psi/dtor c -- profile mixing fraction if (abs(lat).le.15.) then prfrac = 0. latind(1) = 1 latind(2) = 2 else if (abs(lat).le.45.) then prfrac = (abs(lat) - 15.) / 30. latind(1) = 1 latind(2) = 2 else if (abs(lat).le.75.) then prfrac = (abs(lat) - 45.) / 30. latind(1) = 2 latind(2) = 3 else prfrac = 1. latind(1) = 2 latind(2) = 3 endif c -- lambertian surface reflectivity calculation call reflec(iofset,cofs,phi,p1,pr,lat,pteran,rad,ref) c -- flux calculation call flux(indsol,latind,sza,satza,phi,psi,pteran,ozone,ref, 1 rad,minref,res,algflg,snow,surf,cluvf,uvf,eryflux,eryclear) write (6,'(i5,i4,i6,i3,2f7.2,f8.2,f6.3,f6.1,12f8.2,i3)') 1 iyr,iday,itime,isamp,satza,lat,long,r2,ozone,100.*minref, 1 100.*ref,cluvf,eryclear,uvf,eryflux,algflg enddo 900 continue close (24) end