! module co2cooling ! ! Co2 cooling parameterization. Uses V.I. Fomichev May, 1996 ! 15 um co2 band cooling. ! ! subroutine co2cool ! call colum ! (column.F) ! call cco2gr ! Calculate parameterization coeficients for 15 um CO2 band cooling ! call a18int, a18lin (util.F) ! call recur ! call pcoolg ! end subroutine co2cool ! block data pco2o3 ! [init common blocks] ! end block data pco2o3 ! ! Common blocks below are shared with block data pco2o3, which is ! in this file below, but is not part of this module. ! implicit none ! real :: amat,bmat,al common /co2cfg/ amat(43,9),bmat(43,9),al(17) ! ! data from BLOCK DATA PCO2O3 real :: xr,xrh,xrf common /pirgrd/ xr(67), xrh(59), xrf(17) real :: co2vp,co2cl common /pirco2/ co2vp(67), co2cl(67) real :: a150,b150,a360,b360,a540,b540,a720,b720,co2o common /pirlte/ a150(43,9),b150(43,9), a360(43,9),b360(43,9), | a540(43,9),b540(43,9), a720(43,9),b720(43,9), | co2o(4) integer :: ig common /pirmgr/ IG(9) real :: uco2ro,alo,cor150,cor360,cor540,cor720,uco2co common /pirne/ uco2ro(51),alo(51), | cor150(6),cor360(6),cor540(6),cor720(6),uco2co(6) real :: ao3 common /piro3/ ao3(35,9) contains !----------------------------------------------------------------------- ! subroutine co2cool(xlco2e,xlco2d) use params,only: nlev use restart,only: atm use cons,only: co2mix use flds_atmos,only: tn,xno,xno2,xnn2,xnco2,xno3 implicit none ! ! Args: real,intent(out) :: xlco2e(nlev),xlco2d(nlev) ! ! Local: integer :: i,k,kk integer,save :: ncalls=0 real,dimension(59) :: h real,dimension(63) :: hn,hnn real,dimension(67) :: t,x,o3,pt,px,tnn,xmun2,xmuo2,xmuo,xmuo3, | xmuco2,xmuam,xcp,xmden,co2col real,dimension(17) :: co2,o2,sn2,o,am,uco2,pco2,po2,psn2,po,pam, | pden real,dimension(nlev) :: uno,uno2,unn2,unco2,uno3,uam,cpp,uden, | xxn,col,xp real :: ur,xmt,ppn2,ppo,ppo2,cf,dum,xmlw,utop,rco2,zo,flux,cp, | z,zo2,zn2,h1,s,alambda,ak,boltz,a10,y,tt,const ! ncalls = ncalls+1 ! write(6,"('Enter co2cool..')") ! ! Lower atmosphere to zp =-17. is the same as fomichev data ! do i=1,18 ! was do 6 ! ! atm(natm=71,nfatm=8) was read from source history (module restart) x(i) = atm(i,1) t(i) = atm(i,2) o3(i) = atm(i,7) px(i) = x(i) pt(i) = t(i) ! write(6,"('co2cool: i=',i3,' x,t,o3,px,pt=',5e12.4)") ! | i,x(i),t(i),o3(i),px(i),pt(i) enddo ! ! Interpolation between grids ! ur = 8.3144e+7 do k=1,nlev ! was do 7 xmt = xnn2(k)+xno2(k)+xno(k) uno(k) = xno(k)/xmt uno2(k) = xno2(k)/xmt unn2(k) = xnn2(k)/xmt unco2(k) = xnco2(k)/xmt uno3(k) = xno3(k)/xmt uam(k) = (28.*xnn2(k)+32.*xno2(k)+16.*xno(k))/xmt uden(k) = xmt ppn2 = xnn2(k)/xmt ppo = xno(k)/xmt ppo2 = 1.-ppn2-ppo cf = 0.5*ur*(ppo2*7./32.+ppn2*7./28.+ppo*5./16.) cpp(k) = 86400./cf xxn(k) = xnco2(k) enddo ! call fprint8('co2cool',nlev, ! | uno , uno2 , unn2, unco2 ,uno3 , uam , dum,dum, ! | 'uno','uno2','unn2','unco2','uno3','uam',' ' ,' ') xmlw = 44. ! ! Sub colum (column.F): xmlw and xxn(nlev) are input, ! col(nlev) is output. call colum(xmlw,xxn,col) ! write(6,"('co2cool after colum: col=',/,(6e12.4))") col utop = col(49) do k=1,49 ! was do 8 kk = 18+k tnn(kk) = tn(k) xmun2(kk) = unn2(k) xmuo2(kk) = uno2(k) xmuo3(kk) = uno3(k) xmuo(kk) = uno(k) xmuco2(kk)= unco2(k) co2col(kk)= col(k) xmuam(kk) = uam(k) xmden(kk) = uden(k) xcp(kk) = cpp(k) enddo ! write(6,"('co2cool: tnn(19:67) =',/,(6e12.4))") tnn(19:67) ! write(6,"('co2cool: xmuo(19:67) =',/,(6e12.4))") xmuo(19:67) ! write(6,"('co2cool: xmuo2(19:67) =',/,(6e12.4))") xmuo2(19:67) ! write(6,"('co2cool: xmun2(19:67) =',/,(6e12.4))") xmun2(19:67) ! write(6,"('co2cool: xmuco2(19:67)=',/,(6e12.4))") xmuco2(19:67) ! write(6,"('co2cool: xmuo3(19:67) =',/,(6e12.4))") xmuo3(19:67) ! write(6,"('co2cool: xmuam(19:67) =',/,(6e12.4))") xmuam(19:67) ! write(6,"('co2cool: xcp(19:67) =',/,(6e12.4))") xcp(19:67) do k=19,67 ! was do 10 x (k) = atm(k,1) t (k) = tnn(k) o3(k) = xmuo3(k) px(k) = x(k) pt(k) = t(k) ! write(6,"('co2cool: k=',i3,' x,t,o3,px,pt=',5e12.4)") ! | k,x(k),t(k),o3(k),px(k),pt(k) enddo do k=51,67 ! was do 12 kk=k-50 o (kk) = xmuo (k) o2 (kk) = xmuo2 (k) sn2 (kk) = xmun2 (k) co2 (kk) = xmuco2(k) uco2(kk) = co2col(k) am (kk) = xmuam (k) pden(kk) = xmden (k) pco2(kk) = co2 (kk) po2 (kk) = o2 (kk) po (kk) = o (kk) pam (kk) = am (kk) psn2(kk) = sn2 (kk) ! write(6,"('co2cool: kk=',i3,' o,o2,sn2,co2,uco2=',5e12.4)") ! | kk,o(kk),o2(kk),sn2(kk),co2(kk),uco2(kk) enddo ! kk=51,67 do k=51,67 kk=k-50 ! write(6,"('co2cool: kk=',i3,' am,pden,pco2,po2,po=',5e12.4)") ! | kk,am(kk),pden(kk),pco2(kk),po2(kk),po(kk) enddo do k=51,67 kk=k-50 ! write(6,"('co2cool: kk=',i3,' pam,psn2=',2e12.4)") ! | kk,pam(kk),psn2(kk) enddo ! ! CO2 concentration below x=12.5 ! rCO2 = 0.7200E-03 ! rCO2 = 0.3500E-03 ! rCO2 = 0.1500E-03 rco2 = co2mix ! ! To determine the coefficients in both LTE and NLTE regions for a certain ! CO2 concentration. Should be called if CO2 concentration within x=0-16.5 ! has been changed. ! if (ncalls==1) call cco2gr(rco2) ! ! Calculate cooling rate. call recur(uco2) ! zo = 1.5e-12 ! zo = 3.e-12 ! zo = 6.e-12 call pcoolg(h,px,pt,O3,pco2,po2,psn2,po,pam,pden,zo,flux) ! do i=1,59 ! was do 15 ! ! cp - specific heat at constant pressure if (i <= 42) then cp = 7./2.*8.31441e7/28.96 else cp = 8.31441e7/pam(i-42)*(7./2.*(1.-po(i-42))+5./2.*po(i-42)) endif ! ! to convert heating rate (at x=2-12.5) into "K/day" hn(i) = h(i)*86400./cp hnn(i) = h(i) enddo ! i=1,59 ! ! Print results: do k=1,49 ! was do 50 kk = 11+(k-1) xp(k) = x(kk+8) xlco2d(k) = hn(kk) xlco2e(k) = hnn(kk) ! write(6,"('co2cool: k=',i3,' xp,xlco2d,xlco2e=',3e12.4)") ! | k,xp(k),xlco2d(k),xlco2e(k) enddo ! k=1,49 ! do k=50,nlev ! was do 51 cp = 8.31441e7/uam(k)*(7./2.*(1.-uno(k))+5./2.*uno(k)) ! cp = cpp(k) ! Some constants from pcoolg a10 = 1.5988 boltz = 1.38066e-16 ak = 960.217 const = 2.63187e11 ! ! Determine quantum survival probability -- alambda: ! ! --- V-T constants for O2 amd N2 tt=TN(K) y=tt**(-1./3.) zn2=5.5e-17*sqrt(tt)+6.7e-10*exp(-83.8*y) zo2=1.e-15*exp(23.37-230.9*y+564.*y*y) zo=(2.6-0.004*tt)*1.e-12 if (tt < 260.) zo=1.56e-12 if (tt > 260.) zo=1.40e-12 ! ! --- collisional deactivation rate: z= xnn2(k)*zn2+xno2(k)*zo2+xno(k)*zo alambda = a10/(a10+z) ! ! source function s = exp(-ak/TN(K)) ! ! cooling rate h1 = const/uam(k)*unco2(k)*(1.-alambda)*(flux-s) xlco2d(k) = h1*86400./cp xlco2e(k) = h1 ! write(6,"('end co2cool: k=',i3,' xlco2d,xlco2e=',2e12.4)") ! | k,xlco2d(k),xlco2e(k) enddo ! k=50,nlev ! call fprint8('end co2cool',nlev, ! | xp , xlco2d , xlco2e , dum , dum , dum, dum ,dum, ! | 'xp','xlco2d','xlco2e',' ',' ',' ',' ' ,' ') end subroutine co2cool ! C**************************************************************************** C* * C* SUBROUTINE CCO2GR C* * C**************************************************************************** C C to calculate parameterization coeficients for 15 um CO2 band cooling C rate calculation for arbitrary CO2 volume mixing ratio (without height C variation up to x=12.5) in a range of 150-720ppm. C Coefficients are calculated for atmosphere layer from the pressure scale C height x=2 up to x=16.5. Irregular vertical grid is used to account for C internal heat exchange within "LTE layer" (x=0-12.5) C ***** Important!!!!!***** C From x=2 up to x=13.75 coeficients can be calculated only for a vertical C grid with a step of 0.25. Eventually, any vertical step could be used C within the atmospheric layer of x=14-16.5. However, to do it, some C modification is needed. No coeficients are needed to calculate cooling C rates above x=16.5 (see driving programm). C C May, 1996. V.I. Fomichev. C Called by a driving programm C Calls A18LIN (linear interpolation) C Calls A18INT (third order spline) subroutine cco2gr(rco2) ! ! Args: real,intent(in) :: rco2 C input: C RCO2 - CO2 volume mixing in the region below x=12.5 C initial data from BLOCK DATA PCO2O3 come through common blocks C output: parameterization coefficients for both, the matrix parameterization C (is used between x=2 and 12.5) and reccurence formula. C Passing on to other subroutines through common block CO2CFG. C AMAT,BMAT(43,9) - coefficients for the matrix parameterization ! ! Local: integer :: i,j,isgn real :: a real uref(4), co2int(4), uco2(17) ! ! External (util.F): real,external :: a18lin ! ! Calculate coefficients for the matrix paramerization: ! ! subroutine a18int(x1,y1,x2,y2,n1,n2) ! integer,intent(in) :: n1,n2 ! real,intent(in) :: x1(n1),y1(n1) ! real,intent(out) :: x2(n2),y2(n2) ! real function a18lin(x,xn,yn,m,n) ! integer,intent(in) :: m,n ! real,intent(in) :: x,xn(n),yn(n) do 1 i = 1,43 do 1 j = 1,9 if((i.le.5).and.(j.eq.2)) goto 1 isgn = int(sign(1.,a150(i,j))+sign(1.,a360(i,j))+ + sign(1.,a540(i,j))+sign(1.,a720(i,j))) co2int(1)=a150(i,j)/co2o(1) co2int(2)=a360(i,j)/co2o(2) co2int(3)=a540(i,j)/co2o(3) co2int(4)=a720(i,j)/co2o(4) if(isgn.eq.-4) then co2int(1) = alog(-co2int(1)) co2int(2) = alog(-co2int(2)) co2int(3) = alog(-co2int(3)) co2int(4) = alog(-co2int(4)) a = -exp(a18lin(rco2,co2o,co2int,1,4)) else if (isgn.eq.4) then co2int(1) = alog(co2int(1)) co2int(2) = alog(co2int(2)) co2int(3) = alog(co2int(3)) co2int(4) = alog(co2int(4)) a = exp(a18lin(RCO2,co2o,co2int,1,4)) else call a18int(co2o,co2int,RCO2,a,4,1) end if AMAT(i,j)=a*RCO2 isgn = int(sign(1.,b150(i,j))+sign(1.,b360(i,j))+ + sign(1.,b540(i,j))+sign(1.,b720(i,j))) co2int(1)=b150(i,j)/co2o(1) co2int(2)=b360(i,j)/co2o(2) co2int(3)=b540(i,j)/co2o(3) co2int(4)=b720(i,j)/co2o(4) if(isgn.eq.-4) then co2int(1) = alog(-co2int(1)) co2int(2) = alog(-co2int(2)) co2int(3) = alog(-co2int(3)) co2int(4) = alog(-co2int(4)) a = -exp(a18lin(RCO2,co2o,co2int,1,4)) else if (isgn.eq.4) then co2int(1) = alog(co2int(1)) co2int(2) = alog(co2int(2)) co2int(3) = alog(co2int(3)) co2int(4) = alog(co2int(4)) a = exp(a18lin(RCO2,co2o,co2int,1,4)) else call a18int(co2o,co2int,RCO2,a,4,1) end if BMAT(i,j)=a*RCO2 1 continue return end subroutine cco2gr !----------------------------------------------------------------------- SUBROUTINE RECUR(UCO2) C CO2(17) - CO2 vmr at a grid of x=12.5-16.5 with a step = 0.25 C AM(17) - air molecular weight (g/mole) at the same grid C UTOP - CO2 column amount (cm-2) above the level of x=16.5 C AL(17) - coeficients for the reccurence formula. Note that starting up C from x=14.00, these coefficients are identical to escape functions C which can be calculated at any arbitrary vertical grid. C So that starting up from this level any arbitrary vertical grid C can be used. C initial data from BLOCK DATA PCO2O3 come through common blocks ! ! Local: real uref(4), co2int(4), uco2(17) real :: a,cor integer :: i ! ! External (util.F): real,external :: a18lin c c c uco2(17) (CO2 column amount) for CO2. c c calculate coeficients for the reccurence formula: c between x=12.5 and 13.75 these coefficients (al) are calculated using c correction to escape function. Starting up from x=14.00 parameterization c coeficients equal escape function. do 3 i=1,6 a = a18lin(uco2(i),uco2ro,alo,1,51) co2int(1)=cor150(i) co2int(2)=cor360(i) co2int(3)=cor540(i) co2int(4)=cor720(i) uref(1) =uco2co(i)*150./360. uref(2) =uco2co(i) uref(3) =uco2co(i)*540./360. uref(4) =uco2co(i)*720./360. cor = a18lin(uco2(i),uref,co2int,1,4) AL(i)=exp(cor+a) 3 continue do 4 i=7,17 a = a18lin(uco2(i),uco2ro,alo,1,51) AL(i)=exp(a) 4 continue return end subroutine recur C**************************************************************************** C* * C* SUBROUTINE PCOOLG * C* * C**************************************************************************** C C to calculate the heating rate in both, 15 um CO2 and 9.6 um O3 bands. C NLTE conditions are taken into account for CO2 band. Irregular vertical C grid to account for heat exchange in the "LTE-layer" (x=2-12.5), is used. C The subroutine can be used for an arbitrary CO2 volume mixing ratio C with the parameterization coefficients calculated by the CCO2GR C C May, 1996. V.I. Fomichev. C Called by a driving programm C Calls nothing SUBROUTINE PCOOLG(H,X,T,O3,CO2,O2,SN2,O,AM,DEN,ZCO2O,FLUX) C input: C T(67) - temperature (K) at x(=0-16.5, 0.25) levels C X(67) - pressure scale height array; X from 0 up to 16.5, step 0.25 C O3(67) - ozone volume mixing ratio (vmr) at x=0-11., 0.25 C CO2,O2,N2,O(17) - volume mixing ratios (vmr's) for corresponding gases C at X = 12.5 - 16.5, with a step of 0.25 C AM(17) - air molecular weight at X = 12.5 - 16.5, 0.25 C ZCO2O - CO2-O collisional deactivation rate constant (in cm3/sec) C parameterization coeficients for CO2 (calculated by CCO2GR) come from C the CO2CFG and PIRO3 common blocks. Irregular vertical grid (ig(9)) c comes from the BLOCK DATA PCO2O3 (PIRMGR common block) C output: C H(59) - heating rates in erg/g/sec at X = 2-16.5, 0.25 C FLUX - upward flux at X=16.5 to calculate the heating rate above this level. C C Note: 1) Eventually, any arbitrary vertical grid could be utilized starting C up from x=13.75. To do so, the coefficients for the reccurence C formula must be determined at a proper grid (using CCO2GR) and all C input information should be given at this grid too. C 2) As of now, ZCO2O is recommended to vary in a range of (1.5-6)e-12 C with a mean value of 3e-12 cm3/sec, without T-dependence C 3) To get the heating rate values in K/sec, the heating rates C calculated here should be divided by Cp - specific heat at constant C pressure. Where Cp = R/AM * {7/2*SUM[Rv,i - for 2-atom gases] + C 5/2*SUM[Rv,i - for 1-atom gases]}, R - gas constant, Rv - vmr's. C real,dimension(67),intent(in) :: x,t,o3 real,dimension(17),intent(in) :: co2,o2,sn2,o,am,den real,intent(in) :: zco2o real,intent(out) :: h(59),flux C INTERNAL arrays: C su(67) - exponential part of the Planck function for CO2 band at X(67) grid C so3(49) - the same for O3 band at x=(0-11,0.25) and 0 for x=11.25-11.5 C lambda(17) - quantum survival probability at the X grid for levels C starting up from x=12.5 (indexes 51-67) ! ! Local: integer :: i,j,is,js,im real su(67), so3(67), lambda(17) real :: zo2,zn2,zoco2,tt,y,aa1,aa2,h1,h2,h3,d1,d2,d,z ! data so3/49*0./ c const - constant used to determine the heating rate: c const = Na*h*c*V*Gv'/Gv*Avv*1.03, where c Na - Avogadro number, h- Planck's constant, c- light speed, c V - frequency of the fundamental vibrational transition for the c main isotope, Gv'/Gv = 2 - ratio of the statistical weights for the c fundamental transition, Av'v - Einstine coefficient for the c fundamental band of the main isotope, 1.03 - correction to account c for others than fundamental bands in the reccurence formula c constb = const/28.96 - to determine a boundary condition for the reccurence c formula at x=12.5 (air molecule weight is supposed to be 28.96 c below x=12.5... could be changed, if neccessary) c boltz - Boltzmann's constant c a10 - Einstine's coefficient for the fundamental band for the main isotope real :: const,constb,boltz,a10 data const/2.63187e11/, constb/9.08795e9/, : boltz/1.38066e-16/, a10/1.5988/ c aku = h*c/k*Vco2 - to determine the Planck function, k - Boltzmann's constant c ako3 - the same as aku, but for O3 real :: aku,ako3 data aku/-960.217/, ako3/-1500./ so3(:) = 0. c to determine su(67), so3(67): do 1 i = 1,67 su(i)=exp(aku/T(i)) if(i.ge.46) goto 1 so3(i)=exp(ako3/T(i)) 1 continue c to determine lambda(17) do 2 i=1,17 c --- CO2-O2 and CO2-N2 V-T constants: tt=T(i+50) y=tt**(-1./3.) zn2=5.5e-17*sqrt(tt)+6.7e-10*exp(-83.8*y) zo2=1.e-15*exp(23.37-230.9*y+564.*y*y) zoco2=(2.6-0.004*tt)*1.e-12 if (tt < 260.) zoco2=1.56e-12 if (tt > 260.) zoco2=1.40e-12 c --- c air number density: d=den(i) c --- collisional deactivation rate: z=(SN2(i)*zn2+O2(i)*zo2+O(i)*zoco2)*d lambda(i) = a10/(a10+z) 2 continue c cooling rate in both O3 and CO2 bands (x=2-10.5, the matrix approach is used) do 3 i=1,5 is = i+8 h2=(amat(i,1)+bmat(i,1)*su(is))*su(1) h3=ao3(i,1)*so3(1) do 4 j=3,9 js=is+ig(j) h2=h2+(amat(i,j)+bmat(i,j)*su(is))*su(js) 4 h3=h3+ao3(i,j)*so3(js) 3 H(i) = h2+h3*O3(is) do 5 i=6,18 is = i+8 h2=(amat(i,1)+bmat(i,1)*su(is))*su(1) h3=ao3(i,1)*so3(1) do 6 j=2,9 js=is+ig(j) h2=h2+(amat(i,j)+bmat(i,j)*su(is))*su(js) 6 h3=h3+ao3(i,j)*so3(js) 5 H(i) = h2+h3*O3(is) do 7 I=19,35 is=i+8 h2=0. h3=0. do 8 j=1,9 js=is+ig(j) h2=h2+(amat(i,j)+bmat(i,j)*su(is))*su(js) 8 h3=h3+ao3(i,j)*so3(js) 7 H(i) = h2+h3*O3(is) c cooling rate in CO2 bands (x=10.75-12.5, the matrix approach is used) do 9 I=36,43 is=i+8 h2=0. do 10 j=1,9 js=is+ig(j) 10 h2=h2+(amat(i,j)+bmat(i,j)*su(is))*su(js) 9 H(i) = h2 c calculate the heating rate for x=12.75-16.5 (the reccurence formula is used) c --- to form the boundary condition at x=12.5 h1=H(43)/(CO2(1)*(1.-lambda(1))*constb) c --- the reccurence formula do 11 I=2,17 im=i-1 aa1=1.-lambda(im)*(1.-.25*AL(i)-.75*AL(im)) aa2=1.-lambda(i)*(1.-.75*AL(i)-.25*AL(im)) d1=-.25*(AL(i)+3.*AL(im)) d2=.25*(3.*AL(i)+AL(im)) h2=(aa1*h1-d1*su(im+50)-d2*su(i+50))/aa2 H(i+42)=h2*CO2(i)*(1.-lambda(i))/AM(i)*const h1=h2 11 continue c to determine FLUX c cooling rate above x=16.5 is suggested to be calculated by the formula c H(i) = const/AM(i)*CO2(i)*(1.-lambda(i))*(FLUX-su(i)) c no any parameterization coefficients are needed and an arbitrary hight grid c could be used above x=16.5 level FLUX = h2 + su(67) return end subroutine pcoolg end module co2cooling C**************************************************************************** C* * C* BLOCK DATA PCO2O3 * C* * C**************************************************************************** C C to determine the parameterization coefficients for IR CO2 and O3 C parameterizations. Irregular grid for height intergation is used. C C May, 1996. V.I. Fomichev C Called by C Calls nothing BLOCK DATA PCO2O3 c INDEX to organize the first entry calculations common /PIRFEC/ IREP c VERTICAL GRIDs to be used in IR scheme c XR(67) - pressure scale heights, psh's, (=0-16.5) at which input parameter c should be given c XRH(59) - psh's at which cooling rates compute (=2-16.5) c XRF(17) - psh's at which the reccurence formula is utilized (=12.5-16.5) common /PIRGRD/ XR(67),XRH(59),XRF(17) c DATA on CO2 abundance c CO2VP(67) - volume mixing ratio (vmr) for basic CO2 prifile at XR grid c CO2CL(67) - CO2 column amount at XR grid calcilated with use of CO2VP common /PIRCO2/ CO2VP(67), CO2CL(67) c DATA for "LTE" parameterization (x=2-12.5) c IG(9) - indexes of level which should be used to account for the internal c atmospheric heat exchange c AO3(35,9) - coefficients for O3 scheme to calculate cooling rate in c "erg/g/s" at levels x=2-10.5 (with a step of 0.25) c CO2O(4) - vmr for basic CO2 c "LTE-coefficients" for CO2 scheme using to calculate cooling rate in c "erg/g/s" in region x=2-12.5 (with a step of 0.25). To account for internal c heat exchange 9 level in atmosphere are needed. c A150, B150(43,9) - for 150ppm of CO2 c A360, B360(43,9) - for 360ppm of CO2 c A540, B540(43,9) - for 540ppm of CO2 c A720, B720(43,9) - for 720ppm of CO2 common /PIRLTE/ A150(43,9),B150(43,9), A360(43,9),B360(43,9), 1 A540(43,9),B540(43,9), A720(43,9),B720(43,9), 2 CO2O(4) common /PIRMGR/ IG(9) common /PIRO3/ AO3(35,9) c DATA for NLTE parameterization CO2 (x=12.5-16.5) c UCO2RO, ALO(51) - CO2 column amount and corresponding escape functions c (eventually, their log) c COR150, COR360, COR540, COR720(6) - correction to escape functions to c calculate coefficients for the reccurence formula between x=12.5 and 13.75 c UCO2CO(6) - CO2 column amount at x=12.5-13.75 (step - 0.25) for 360 ppm c to be used to correct escape functions in this region common /PIRNE/ UCO2RO(51),ALO(51), 1 COR150(6),COR360(6),COR540(6),COR720(6),UCO2CO(6) data IREP/0/ data (XR(i), i=1,67)/ & 0.0, 0.25, 0.5, 0.75, 1.0, 1.25, 1.5, 1.75, 2.0, 2.25, 2.5, 2.75, & 3.0, 3.25, 3.5, 3.75, 4.0, 4.25, 4.5, 4.75, 5.0, 5.25, 5.5, 5.75, & 6.0, 6.25, 6.5, 6.75, 7.0, 7.25, 7.5, 7.75, 8.0, 8.25, 8.5, 8.75, & 9.0, 9.25, 9.5, 9.75,10.0,10.25,10.5,10.75,11.0,11.25,11.5,11.75, &12.0,12.25,12.5,12.75,13.0,13.25,13.5,13.75,14.0,14.25,14.5,14.75, &15.0,15.25,15.5,15.75,16.0,16.25,16.5/ data (XRH(i), i=1,59)/ & 2.0, 2.25, 2.5, 2.75, 3.0, 3.25, 3.5, 3.75, 4.0, 4.25, 4.5, 4.75, & 5.0, 5.25, 5.5, 5.75, 6.0, 6.25, 6.5, 6.75, 7.0, 7.25, 7.5, 7.75, & 8.0, 8.25, 8.5, 8.75, 9.0, 9.25, 9.5, 9.75,10.0,10.25,10.5,10.75, &11.0,11.25,11.5,11.75,12.0,12.25,12.5,12.75,13.0,13.25,13.5,13.75, &14.0,14.25,14.5,14.75,15.0,15.25,15.5,15.75,16.0,16.25,16.5/ data (XRF(i), i=1,17)/ &12.5,12.75,13.0,13.25,13.5,13.75,14.0,14.25,14.5,14.75,15.0,15.25, &15.5,15.75,16.0,16.25,16.5/ data (CO2VP(i), i=1,67) / 50*3.600E-04, & 3.58E-04, 3.54E-04, 3.50E-04, 3.41E-04, 3.28E-04, 3.11E-04, & 2.93E-04, 2.75E-04, 2.56E-04, 2.37E-04, 2.18E-04, 1.99E-04, & 1.80E-04, 1.61E-04, 1.42E-04, 1.24E-04, 1.06E-04/ data (CO2CL(i), i=1,67)/ & 7.760162E+21,6.043619E+21,4.706775E+21,3.665639E+21,2.854802E+21, & 2.223321E+21,1.731524E+21,1.348511E+21,1.050221E+21,8.179120E+20, & 6.369897E+20,4.960873E+20,3.863524E+20,3.008908E+20,2.343332E+20, & 1.824981E+20,1.421289E+20,1.106894E+20,8.620418E+19,6.713512E+19, & 5.228412E+19,4.071814E+19,3.171055E+19,2.469544E+19,1.923206E+19, & 1.497717E+19,1.166347E+19,9.082750E+18,7.072886E+18,5.507602E+18, & 4.288557E+18,3.339164E+18,2.599777E+18,2.023941E+18,1.575480E+18, & 1.226217E+18,9.542118E+17,7.423736E+17,5.773939E+17,4.489076E+17, & 3.488423E+17,2.709114E+17,2.102188E+17,1.629513E+17,1.261393E+17, & 9.747013E+16,7.514254E+16,5.775381E+16,4.421144E+16,3.366464E+16, & 2.546953E+16,1.913609E+16,1.425730E+16,1.052205E+16,7.700499E+15, & 5.596946E+15,4.044901E+15,2.905406E+15,2.073209E+15,1.469497E+15, & 1.033831E+15,7.212717E+14,4.986668E+14,3.415852E+14,2.320320E+14, & 1.565317E+14,1.052353E+14/ data (IG(i),i=1,9)/-25,-12,-7,-3,-1,0,1,3,6/ data (CO2O(i),i=1,4)/150.e-6, 360.e-6, 540.e-6, 720.e-6/ data (UCO2CO(i),i=1,6) /2.546953E+16,1.913609E+16,1.425730E+16, & 1.052205E+16,7.700499E+15,5.596946E+15/ data (COR150(i),i=1,6) /2.530798E-01,2.048590E-01,1.050616E-01, & 9.172653E-02,4.528593E-02,3.384089E-02/ data (COR360(i),i=1,6) /4.848122E-01,4.224252E-01,2.958043E-01, & 2.067767E-01,1.244037E-01,5.792163E-02/ data (COR540(i),i=1,6) /6.263170E-01,5.332695E-01,3.773434E-01, & 2.592216E-01,1.514242E-01,6.889337E-02/ data (COR720(i),i=1,6) /7.248522E-01,6.199704E-01,4.525580E-01, & 3.046277E-01,1.787235E-01,7.953181E-02/ data (UCO2RO(i),i=1,51) /2.699726E+11,5.810773E+11,1.106722E+12, & 1.952319E+12,3.306797E+12,5.480155E+12,8.858565E+12,1.390142E+13, & 2.129301E+13,3.209300E+13,4.784654E+13,7.091442E+13,1.052353E+14, & 1.565317E+14,2.320320E+14,3.415852E+14,4.986668E+14,7.212717E+14, & 1.033831E+15,1.469497E+15,2.073209E+15,2.905406E+15,4.044901E+15, & 5.596946E+15,7.700499E+15,1.052205E+16,1.425730E+16,1.913609E+16, & 2.546953E+16,3.366464E+16,4.421144E+16,5.775381E+16,7.514254E+16, & 9.747013E+16,1.261393E+17,1.629513E+17,2.102188E+17,2.709114E+17, & 3.488423E+17,4.489076E+17,5.773939E+17,7.423736E+17,9.542118E+17, & 1.226217E+18,1.575480E+18,2.023941E+18,2.599777E+18,3.339164E+18, & 4.288557E+18,5.507602E+18,7.072886E+18/ data (ALO(i),i=1,51) /-2.410106E-04,-5.471415E-04,-1.061586E-03, & -1.879789E-03,-3.166020E-03,-5.185436E-03,-8.216667E-03, & -1.250894E-02,-1.838597E-02,-2.631114E-02,-3.688185E-02, & -5.096491E-02,-7.004056E-02,-9.603746E-02,-1.307683E-01, & -1.762946E-01,-2.350226E-01,-3.095215E-01,-4.027339E-01, & -5.178570E-01,-6.581256E-01,-8.265003E-01,-1.024684E+00, & -1.252904E+00,-1.509470E+00,-1.788571E+00,-2.081700E+00, & -2.379480E+00,-2.675720E+00,-2.967325E+00,-3.252122E+00, & -3.530485E+00,-3.803720E+00,-4.072755E+00,-4.338308E+00, & -4.601048E+00,-4.861585E+00,-5.120370E+00,-5.377789E+00, & -5.634115E+00,-5.889388E+00,-6.143488E+00,-6.396436E+00, & -6.648774E+00,-6.901465E+00,-7.155207E+00,-7.409651E+00, & -7.663536E+00,-7.915682E+00,-8.165871E+00,-8.415016E+00/ data ((AO3(ix,js),js=1,9),ix=1,10)/ & 5.690E+09, 0.000E+00, 1.962E+09, 5.015E+09, 7.105E+09,-3.952E+10, & 8.916E+09, 1.204E+09, 5.434E+09, 4.997E+09, 0.000E+00, 2.168E+09, & 4.981E+09, 7.327E+09,-3.855E+10, 9.296E+09, 1.326E+09, 4.575E+09, & 4.179E+09, 0.000E+00, 2.198E+09, 4.594E+09, 7.780E+09,-3.689E+10, & 9.624E+09, 1.107E+09, 3.736E+09, 3.469E+09, 0.000E+00, 1.941E+09, & 3.921E+09, 8.423E+09,-3.508E+10, 9.859E+09, 8.736E+08, 3.058E+09, & 2.930E+09, 0.000E+00, 1.649E+09, 3.378E+09, 8.985E+09,-3.362E+10, & 9.941E+09, 8.110E+08, 2.367E+09, 2.020E+09, 3.866E+08, 1.368E+09, & 3.049E+09, 9.452E+09,-3.249E+10, 1.010E+10, 6.212E+08, 2.117E+09, & 1.652E+09, 3.711E+08, 1.175E+09, 2.747E+09, 9.818E+09,-3.165E+10, & 1.014E+10, 5.470E+08, 1.843E+09, 1.391E+09, 3.497E+08, 1.033E+09, & 2.570E+09, 1.008E+10,-3.102E+10, 1.017E+10, 5.121E+08, 1.552E+09, & 1.197E+09, 3.329E+08, 9.689E+08, 2.427E+09, 1.024E+10,-3.056E+10, & 1.011E+10, 4.521E+08, 1.312E+09, 9.961E+08, 3.390E+08, 9.867E+08, & 2.375E+09, 1.032E+10,-3.033E+10, 1.011E+10, 5.322E+08, 1.178E+09/ data ((AO3(ix,js),js=1,9),ix=11,20)/ & 8.408E+08, 4.922E+08, 9.226E+08, 2.358E+09, 1.039E+10,-3.030E+10, & 1.004E+10, 5.551E+08, 1.141E+09, 7.250E+08, 5.353E+08, 9.493E+08, & 2.497E+09, 1.040E+10,-3.047E+10, 9.971E+09, 8.389E+08, 9.122E+08, & 6.484E+08, 5.214E+08, 1.014E+09, 2.625E+09, 1.038E+10,-3.058E+10, & 9.903E+09, 7.720E+08, 9.455E+08, 5.785E+08, 4.987E+08, 1.037E+09, & 2.694E+09, 1.039E+10,-3.089E+10, 9.698E+09, 1.068E+09, 8.903E+08, & 5.272E+08, 4.955E+08, 1.194E+09, 2.998E+09, 1.026E+10,-3.160E+10, & 9.466E+09, 1.239E+09, 9.659E+08, 4.891E+08, 5.131E+08, 1.299E+09, & 3.456E+09, 1.016E+10,-3.275E+10, 9.003E+09, 1.593E+09, 1.083E+09, & 4.643E+08, 5.668E+08, 1.500E+09, 4.096E+09, 9.972E+09,-3.422E+10, & 8.345E+09, 1.790E+09, 1.256E+09, 4.461E+08, 6.391E+08, 1.870E+09, & 4.881E+09, 9.608E+09,-3.621E+10, 7.386E+09, 2.105E+09, 1.460E+09, & 6.471E+08, 1.536E+08, 2.567E+09, 5.909E+09, 9.088E+09,-3.844E+10, & 6.248E+09, 2.173E+09, 1.536E+09, 5.377E+08, 1.101E+08, 3.474E+09, & 6.992E+09, 8.333E+09,-4.065E+10, 4.979E+09, 2.019E+09, 1.594E+09/ data ((AO3(ix,js),js=1,9),ix=21,30)/ & 7.064E+08, 1.346E+08, 4.535E+09, 7.814E+09, 7.361E+09,-4.248E+10, & 3.786E+09, 1.733E+09, 1.431E+09, 1.210E+09, 2.133E+08, 5.772E+09, & 8.402E+09, 6.184E+09,-4.400E+10, 2.883E+09, 1.530E+09, 1.163E+09, & 2.121E+09, 3.450E+08, 7.111E+09, 8.611E+09, 4.874E+09,-4.499E+10, & 2.055E+09, 1.126E+09, 8.840E+08, 3.393E+09, 5.066E+08, 8.176E+09, & 8.419E+09, 3.777E+09,-4.573E+10, 1.418E+09, 8.927E+08, 6.356E+08, & 4.838E+09, 6.534E+08, 9.188E+09, 7.895E+09, 2.746E+09,-4.628E+10, & 9.548E+08, 6.033E+08, 4.810E+08, 5.852E+09, 1.286E+09, 9.867E+09, & 7.121E+09, 2.046E+09,-4.663E+10, 6.643E+08, 4.034E+08, 2.258E+08, & 6.624E+09, 2.439E+09, 1.044E+10, 5.950E+09, 1.380E+09,-4.691E+10, & 4.209E+08, 2.635E+08, 2.211E+08, 7.030E+09, 3.865E+09, 1.055E+10, & 4.781E+09, 1.009E+09,-4.708E+10, 3.238E+08, 2.107E+08, 1.094E+08, & 8.099E+09, 4.843E+09, 1.065E+10, 3.732E+09, 6.798E+08,-4.719E+10, & 1.852E+08, 1.233E+08, 1.513E+08, 8.382E+09, 6.426E+09, 1.018E+10, & 2.771E+09, 4.593E+08,-4.726E+10, 1.084E+08, 8.368E+07, 0.000E+00/ data ((AO3(ix,js),js=1,9),ix=31,35)/ & 9.545E+09, 7.667E+09, 9.429E+09, 1.927E+09, 3.143E+08,-4.731E+10, & 7.075E+07, 5.179E+07, 0.000E+00, 1.036E+10, 9.146E+09, 8.014E+09, & 1.355E+09, 2.280E+08,-4.734E+10, 4.459E+07, 3.463E+07, 0.000E+00, & 1.180E+10, 1.021E+10, 6.725E+09, 9.575E+08, 1.091E+08,-4.736E+10, & 3.018E+07, 2.000E+07, 0.000E+00, 8.878E+09, 7.184E+09, 3.569E+09, & 4.023E+08, 4.656E+07,-3.079E+10, 1.374E+07, 0.000E+00, 0.000E+00, & 4.652E+09, 3.469E+09, 1.385E+09, 1.114E+08, 1.368E+07,-1.421E+10, & 3.553E+06, 0.000E+00, 0.000E+00/ data ((A150(ix,js),js=1,9),ix=1,5)/ & 2.5035E+01, 0.0000E+00, 9.4137E+01, 2.4765E+03, 3.2467E+03, &-1.1090E+04, 4.0701E+02, 2.1062E+03, 6.1997E+02, 3.1497E+01, & 0.0000E+00, 1.5661E+02, 3.1886E+03, 3.0556E+03,-1.2560E+04, & 2.1713E+02, 2.2157E+03, 7.7184E+02, 3.6209E+01, 0.0000E+00, & 2.5391E+02, 3.8867E+03, 2.8525E+03,-1.4009E+04,-3.8903E+00, & 2.2576E+03, 9.1027E+02, 3.9456E+01, 0.0000E+00, 4.0539E+02, & 4.4674E+03, 2.6199E+03,-1.5203E+04,-2.0232E+02, 2.2454E+03, & 1.0123E+03, 4.3855E+01, 0.0000E+00, 6.3518E+02, 4.9490E+03, & 2.3873E+03,-1.6282E+04,-4.0148E+02, 2.2334E+03, 1.0983E+03/ data ((A150(ix,js),js=1,9),ix=6,10)/ & 1.4880E+01, 5.2053E+01, 8.8187E+02, 5.4055E+03, 2.1311E+03, &-1.7305E+04,-5.7150E+02, 2.2465E+03, 1.1736E+03, 1.8779E+01, & 8.2377E+01, 1.1180E+03, 5.7943E+03, 1.8908E+03,-1.8247E+04, &-7.1038E+02, 2.3030E+03, 1.2335E+03, 2.3004E+01, 1.2959E+02, & 1.3044E+03, 6.1369E+03, 1.6879E+03,-1.9121E+04,-8.4776E+02, & 2.4019E+03, 1.2764E+03, 2.7051E+01, 1.9749E+02, 1.4427E+03, & 6.4390E+03, 1.4997E+03,-1.9896E+04,-9.8529E+02, 2.5405E+03, & 1.2946E+03, 3.0107E+01, 2.8477E+02, 1.5605E+03, 6.6846E+03, & 1.3469E+03,-2.0582E+04,-1.1115E+03, 2.7010E+03, 1.2951E+03/ data ((A150(ix,js),js=1,9),ix=11,15)/ & 3.3100E+01, 3.8721E+02, 1.6772E+03, 6.8908E+03, 1.2666E+03, &-2.1299E+04,-1.2333E+03, 2.9002E+03, 1.2922E+03, 3.6769E+01, & 5.0026E+02, 1.8186E+03, 7.0844E+03, 1.2377E+03,-2.2143E+04, &-1.3997E+03, 3.1384E+03, 1.3038E+03, 4.0958E+01, 6.2580E+02, & 1.9791E+03, 7.3078E+03, 1.1888E+03,-2.3138E+04,-1.5839E+03, & 3.4194E+03, 1.3343E+03, 4.5506E+01, 7.5175E+02, 2.1209E+03, & 7.5851E+03, 1.1225E+03,-2.4221E+04,-1.8076E+03, 3.7350E+03, & 1.3778E+03, 5.0606E+01, 8.8360E+02, 2.2346E+03, 7.9275E+03, & 1.0078E+03,-2.5290E+04,-2.0653E+03, 4.0577E+03, 1.4264E+03/ data ((A150(ix,js),js=1,9),ix=16,20)/ & 5.6331E+01, 1.0256E+03, 2.3251E+03, 8.3240E+03, 8.5705E+02, &-2.6376E+04,-2.3440E+03, 4.3896E+03, 1.4827E+03, 6.2586E+01, & 1.1799E+03, 2.4192E+03, 8.7537E+03, 7.1111E+02,-2.7532E+04, &-2.6284E+03, 4.7408E+03, 1.5489E+03, 6.9499E+01, 1.3513E+03, & 2.5385E+03, 9.2024E+03, 5.9855E+02,-2.8837E+04,-2.9066E+03, & 5.1230E+03, 1.6315E+03, 9.8718E+01, 1.4643E+03, 2.7198E+03, & 9.6727E+03, 5.6542E+02,-3.0324E+04,-3.1708E+03, 5.5386E+03, & 1.7411E+03, 1.4214E+02, 1.5587E+03, 2.9412E+03, 1.0183E+04, & 6.2959E+02,-3.2074E+04,-3.4238E+03, 5.9832E+03, 1.8906E+03/ data ((A150(ix,js),js=1,9),ix=21,25)/ & 2.0882E+02, 1.6315E+03, 3.1941E+03, 1.0757E+04, 8.0531E+02, &-3.4141E+04,-3.6662E+03, 6.4368E+03, 2.0964E+03, 3.0513E+02, & 1.6987E+03, 3.4635E+03, 1.1429E+04, 1.0747E+03,-3.6585E+04, &-3.8843E+03, 6.8932E+03, 2.3807E+03, 4.3086E+02, 1.7634E+03, & 3.7650E+03, 1.2253E+04, 1.4449E+03,-3.9487E+04,-4.0949E+03, & 7.3521E+03, 2.7604E+03, 5.6228E+02, 1.8545E+03, 4.0922E+03, & 1.3286E+04, 1.8655E+03,-4.2969E+04,-4.3141E+03, 7.7726E+03, & 3.3137E+03, 6.5591E+02, 1.9789E+03, 4.4192E+03, 1.4473E+04, & 2.2732E+03,-4.6758E+04,-4.5119E+03, 8.1060E+03, 4.0069E+03/ data ((A150(ix,js),js=1,9),ix=26,31)/ & 7.1489E+02, 2.1391E+03, 4.7555E+03, 1.5792E+04, 2.6329E+03, &-5.0784E+04,-4.6936E+03, 8.3781E+03, 4.7906E+03, 7.5917E+02, & 2.3306E+03, 5.1321E+03, 1.7236E+04, 2.9779E+03,-5.5179E+04, &-4.8828E+03, 8.7136E+03, 5.5138E+03, 8.1061E+02, 2.5462E+03, & 5.5725E+03, 1.8761E+04, 3.3368E+03,-6.0007E+04,-5.0788E+03, & 9.1950E+03, 6.1026E+03, 8.6939E+02, 2.7684E+03, 6.1037E+03, & 2.0165E+04, 3.7711E+03,-6.4792E+04,-5.2866E+03, 9.9049E+03, & 6.1727E+03, 9.4294E+02, 3.0123E+03, 6.6875E+03, 2.1225E+04, & 4.4266E+03,-6.9341E+04,-5.5392E+03, 1.1040E+04, 5.5220E+03, & 1.0189E+03, 3.2497E+03, 7.2300E+03, 2.0698E+04, 5.1305E+03, &-7.1191E+04,-5.7478E+03, 1.2510E+04, 3.5534E+03/ data ((A150(ix,js),js=1,9),ix=32,37)/ & 1.1449E+03, 3.7084E+03, 7.4630E+03, 1.9495E+04, 5.5511E+03, &-7.3819E+04,-5.4997E+03, 1.4972E+04, 2.3195E+03, 1.2999E+03, & 4.3325E+03, 7.4755E+03, 1.6338E+04, 5.1750E+03,-7.5166E+04, &-4.4753E+03, 1.7530E+04, 2.4555E+03, 1.4512E+03, 5.0331E+03, & 7.0224E+03, 1.5180E+04, 3.2675E+03,-7.8134E+04,-3.0385E+03, & 2.0430E+04, 2.7439E+03, 1.5838E+03, 5.7655E+03, 6.0092E+03, & 1.4199E+04, 5.0935E+03,-8.0938E+04,-3.1873E+03, 2.1397E+04, & 2.3826E+03, 1.6805E+03, 6.6768E+03, 4.2956E+03, 1.4811E+04, & 8.4327E+03,-8.6256E+04,-4.1912E+03, 2.2919E+04, 1.4885E+03, & 1.7701E+03, 7.7803E+03, 3.4348E+03, 1.8569E+04, 1.3687E+04, &-1.0164E+05,-6.1459E+03, 2.7017E+04, 1.4367E+03/ data ((A150(ix,js),js=1,9),ix=38,43)/ & 2.1081E+03, 9.8559E+03, 3.3385E+03, 2.5453E+04, 1.9163E+04, &-1.3508E+05,-7.9100E+03, 4.0371E+04, 3.1637E+03, 2.5117E+03, & 1.2043E+04, 2.9837E+03, 3.0525E+04, 2.1190E+04,-1.6547E+05, &-1.0361E+04, 5.5753E+04, 6.3850E+03, 2.7271E+03, 1.3131E+04, & 1.4005E+03, 3.1777E+04, 1.8261E+04,-1.7548E+05,-1.2891E+04, & 6.4744E+04, 9.9029E+03, 2.4817E+03, 1.2549E+04, 1.2521E+03, & 3.2506E+04, 1.5901E+04,-1.6344E+05,-1.7080E+04, 5.9506E+04, & 7.9804E+03, 2.0450E+03, 1.1403E+04, 3.8851E+03, 3.7098E+04, & 1.3215E+04,-1.5453E+05,-2.1257E+04, 5.2256E+04, 5.4705E+03, & 1.7609E+03, 1.1782E+04, 7.3394E+03, 3.6003E+04, 5.7292E+03, &-1.4489E+05,-2.1518E+04, 4.7419E+04, 2.0765E+03/ data ((B150(ix,js),js=1,9),ix=1,5)/ & 3.5399E+03, 0.0000E+00, 6.0508E+03, 5.9109E+04, 2.4118E+04, &-1.8103E+05,-7.3845E+02, 1.8008E+04, 1.7036E+04, 3.1195E+03, & 0.0000E+00, 6.9773E+03, 5.7201E+04, 2.0239E+04,-1.7100E+05, &-2.3973E+03, 1.8321E+04, 1.5737E+04, 3.0408E+03, 0.0000E+00, & 9.0810E+03, 5.5104E+04, 1.3930E+04,-1.6106E+05,-4.9406E+03, & 1.6082E+04, 1.4928E+04, 3.2990E+03, 0.0000E+00, 1.3161E+04, & 5.6155E+04, 1.0436E+04,-1.6715E+05,-7.0673E+03, 1.5202E+04, & 1.5529E+04, 3.5043E+03, 0.0000E+00, 1.8035E+04, 5.4916E+04, & 7.8353E+03,-1.7025E+05,-8.4919E+03, 1.4317E+04, 1.5764E+04/ data ((B150(ix,js),js=1,9),ix=6,10)/ & 1.8653E+03, 2.7059E+03, 1.9745E+04, 5.4113E+04, 5.9014E+03, &-1.7121E+05,-9.3646E+03, 1.3888E+04, 1.5667E+04, 2.0888E+03, & 3.6216E+03, 1.9181E+04, 5.2430E+04, 4.3250E+03,-1.6730E+05, &-9.8248E+03, 1.3662E+04, 1.5012E+04, 2.2582E+03, 4.7714E+03, & 1.7132E+04, 4.9716E+04, 2.5383E+03,-1.5818E+05,-1.0089E+04, & 1.3299E+04, 1.3889E+04, 2.4496E+03, 6.2997E+03, 1.5321E+04, & 4.8636E+04, 9.3940E+02,-1.5318E+05,-1.0680E+04, 1.3768E+04, & 1.2975E+04, 2.7048E+03, 8.2966E+03, 1.4669E+04, 4.9707E+04, &-6.8884E+02,-1.5464E+05,-1.1834E+04, 1.4872E+04, 1.2636E+04/ data ((B150(ix,js),js=1,9),ix=11,15)/ & 3.0080E+03, 1.0337E+04, 1.4758E+04, 5.1024E+04,-2.2171E+03, &-1.5788E+05,-1.3463E+04, 1.6569E+04, 1.2308E+04, 3.3257E+03, & 1.2025E+04, 1.5122E+04, 5.1259E+04,-3.7550E+03,-1.5938E+05, &-1.5169E+04, 1.8104E+04, 1.1931E+04, 3.6784E+03, 1.3370E+04, & 1.5491E+04, 5.0380E+04,-5.6000E+03,-1.5854E+05,-1.6538E+04, & 1.9311E+04, 1.1575E+04, 4.0171E+03, 1.4413E+04, 1.5479E+04, & 4.8944E+04,-7.1724E+03,-1.5614E+05,-1.7640E+04, 2.0517E+04, & 1.1214E+04, 4.3161E+03, 1.5290E+04, 1.5284E+04, 4.8663E+04, &-7.7645E+03,-1.5691E+05,-1.8632E+04, 2.2492E+04, 1.1169E+04/ data ((B150(ix,js),js=1,9),ix=16,20)/ & 4.6114E+03, 1.6154E+04, 1.4979E+04, 4.9364E+04,-7.9418E+03, &-1.6022E+05,-1.9786E+04, 2.5014E+04, 1.1385E+04, 4.8886E+03, & 1.7121E+04, 1.4681E+04, 5.0603E+04,-8.0295E+03,-1.6487E+05, &-2.1104E+04, 2.7565E+04, 1.1912E+04, 5.1220E+03, 1.8151E+04, & 1.4554E+04, 5.2649E+04,-8.0370E+03,-1.7162E+05,-2.2730E+04, & 3.0243E+04, 1.2705E+04, 6.9450E+03, 1.4807E+04, 1.6724E+04, & 5.6069E+04,-8.0590E+03,-1.8189E+05,-2.4936E+04, 3.3115E+04, & 1.4017E+04, 9.4904E+03, 1.0887E+04, 1.9534E+04, 6.0737E+04, &-8.1973E+03,-1.9502E+05,-2.7754E+04, 3.5773E+04, 1.5904E+04/ data ((B150(ix,js),js=1,9),ix=21,25)/ & 1.3201E+04, 6.7672E+03, 2.2982E+04, 6.6827E+04,-8.7998E+03, &-2.1116E+05,-3.1231E+04, 3.7848E+04, 1.8559E+04, 1.8353E+04, & 2.9582E+03, 2.7018E+04, 7.4757E+04,-1.0238E+04,-2.3087E+05, &-3.5386E+04, 3.9278E+04, 2.2035E+04, 2.4236E+04,-7.1673E+00, & 3.1406E+04, 8.4123E+04,-1.2597E+04,-2.5210E+05,-3.9944E+04, & 3.9900E+04, 2.5639E+04, 2.8899E+04,-1.8050E+03, 3.5598E+04, & 9.3624E+04,-1.6264E+04,-2.7076E+05,-4.3916E+04, 3.8269E+04, & 2.9577E+04, 3.2228E+04,-2.3769E+03, 4.1189E+04, 1.0710E+05, &-2.1940E+04,-2.9820E+05,-4.9133E+04, 3.6912E+04, 3.4890E+04/ data ((B150(ix,js),js=1,9),ix=26,31)/ & 3.4901E+04,-1.7054E+03, 4.9616E+04, 1.2729E+05,-2.9857E+04, &-3.4296E+05,-5.6738E+04, 3.7392E+04, 4.2347E+04, 3.7694E+04, & 2.6968E+02, 6.1729E+04, 1.5470E+05,-4.0354E+04,-4.0768E+05, &-6.7248E+04, 4.0115E+04, 5.1589E+04, 4.0809E+04, 4.1844E+03, & 7.8173E+04, 1.8915E+05,-5.3538E+04,-4.9425E+05,-8.0910E+04, & 4.5420E+04, 6.2215E+04, 4.3379E+04, 1.0525E+04, 9.9171E+04, & 2.2658E+05,-6.7757E+04,-5.9481E+05,-9.6665E+04, 5.4277E+04, & 6.9248E+04, 4.5667E+04, 1.9161E+04, 1.2399E+05, 2.5979E+05, &-8.0920E+04,-6.9722E+05,-1.1391E+05, 6.8928E+04, 6.7840E+04, & 4.8149E+04, 2.9647E+04, 1.5197E+05, 2.6905E+05,-9.0125E+04, &-7.6851E+05,-1.3001E+05, 8.9270E+04, 4.9048E+04/ data ((B150(ix,js),js=1,9),ix=32,37)/ & 5.1966E+04, 4.5702E+04, 1.7720E+05, 2.6058E+05,-9.9247E+04, &-8.4400E+05,-1.4248E+05, 1.1906E+05, 3.5837E+04, 5.7382E+04, & 6.7581E+04, 1.9871E+05, 2.1424E+05,-1.0592E+05,-9.0196E+05, &-1.4324E+05, 1.4835E+05, 3.8639E+04, 6.3804E+04, 9.4177E+04, & 2.0441E+05, 1.9644E+05,-1.3089E+05,-9.7326E+05,-1.3977E+05, & 1.7913E+05, 4.3016E+04, 7.0638E+04, 1.2395E+05, 1.8916E+05, & 1.7195E+05,-1.1044E+05,-1.0265E+06,-1.4544E+05, 1.9545E+05, & 3.6843E+04, 7.6839E+04, 1.6172E+05, 1.4909E+05, 1.6911E+05, &-7.5704E+04,-1.0973E+06,-1.5936E+05, 2.1662E+05, 2.4064E+04, & 8.3232E+04, 2.0819E+05, 1.2835E+05, 2.0677E+05,-3.6307E+04, &-1.2961E+06,-1.9272E+05, 2.6288E+05, 2.1977E+04/ data ((B150(ix,js),js=1,9),ix=38,43)/ & 9.6157E+04, 2.8562E+05, 1.1966E+05, 2.8749E+05,-1.2237E+04, &-1.7528E+06,-2.4263E+05, 4.0970E+05, 4.1317E+04, 1.1357E+05, & 3.7596E+05, 9.8652E+04, 3.6598E+05, 1.6096E+04,-2.2962E+06, &-2.8607E+05, 6.2929E+05, 7.4848E+04, 1.2829E+05, 4.4149E+05, & 4.2542E+04, 4.3449E+05, 3.5764E+04,-2.7327E+06,-3.1740E+05, & 8.5228E+05, 1.2250E+05, 1.2481E+05, 4.3349E+05, 3.3622E+04, & 4.8685E+05, 7.3970E+04,-2.7619E+06,-3.5847E+05, 8.7436E+05, & 1.0368E+05, 1.0994E+05, 3.7480E+05, 1.2328E+05, 6.1220E+05, & 7.4486E+04,-2.8231E+06,-4.3841E+05, 8.4803E+05, 7.9803E+04, & 1.0138E+05, 3.7276E+05, 2.4376E+05, 7.1193E+05, 1.4092E+05, &-3.3243E+06,-4.4325E+05, 1.0532E+06, 3.2305E+04/ data ((A360(ix,js),js=1,9),ix=1,5)/ & 1.7094E+01, 0.0000E+00, 8.2972E+01, 2.5829E+03, 4.4982E+03, &-1.3393E+04, 9.4326E+02, 2.8560E+03, 5.4669E+02, 2.2911E+01, & 0.0000E+00, 1.4184E+02, 3.4315E+03, 4.6065E+03,-1.5742E+04, & 7.9383E+02, 3.2302E+03, 7.5442E+02, 2.7345E+01, 0.0000E+00, & 2.3163E+02, 4.3474E+03, 4.7106E+03,-1.8243E+04, 5.5670E+02, & 3.4688E+03, 9.8490E+02, 3.0814E+01, 0.0000E+00, 3.7251E+02, & 5.2547E+03, 4.7045E+03,-2.0572E+04, 3.1055E+02, 3.5462E+03, & 1.2007E+03, 3.5494E+01, 0.0000E+00, 5.9797E+02, 6.1381E+03, & 4.5752E+03,-2.2700E+04,-9.7490E-02, 3.5437E+03, 1.3996E+03/ data ((A360(ix,js),js=1,9),ix=6,10)/ & 9.1340E+00, 4.8891E+01, 8.6264E+02, 7.0522E+03, 4.2702E+03, &-2.4662E+04,-3.0600E+02, 3.5640E+03, 1.5687E+03, 1.2880E+01, & 8.0570E+01, 1.1535E+03, 7.8722E+03, 3.8761E+03,-2.6396E+04, &-5.6181E+02, 3.6560E+03, 1.7011E+03, 1.7276E+01, 1.3152E+02, & 1.4156E+03, 8.5681E+03, 3.5175E+03,-2.7925E+04,-8.2205E+02, & 3.8262E+03, 1.7943E+03, 2.1698E+01, 2.0535E+02, 1.6357E+03, & 9.1928E+03, 3.1758E+03,-2.9308E+04,-1.0872E+03, 4.0624E+03, & 1.8441E+03, 2.5396E+01, 2.9819E+02, 1.8384E+03, 9.7178E+03, & 2.8742E+03,-3.0523E+04,-1.3368E+03, 4.3248E+03, 1.8525E+03/ data ((A360(ix,js),js=1,9),ix=11,15)/ & 2.8749E+01, 4.0334E+02, 2.0531E+03, 1.0151E+04, 2.6704E+03, &-3.1682E+04,-1.5687E+03, 4.6207E+03, 1.8353E+03, 3.2301E+01, & 5.2404E+02, 2.3304E+03, 1.0546E+04, 2.5451E+03,-3.2994E+04, &-1.8650E+03, 4.9665E+03, 1.8239E+03, 3.6168E+01, 6.6265E+02, & 2.6543E+03, 1.0920E+04, 2.3793E+03,-3.4430E+04,-2.1422E+03, & 5.3802E+03, 1.8279E+03, 4.0172E+01, 8.1883E+02, 2.9511E+03, & 1.1315E+04, 2.2543E+03,-3.6037E+04,-2.4455E+03, 5.8820E+03, & 1.8575E+03, 4.4122E+01, 9.9436E+02, 3.1962E+03, 1.1812E+04, & 2.1571E+03,-3.7838E+04,-2.7947E+03, 6.4488E+03, 1.9130E+03/ data ((A360(ix,js),js=1,9),ix=16,20)/ & 4.8635E+01, 1.1934E+03, 3.3757E+03, 1.2365E+04, 2.0383E+03, &-3.9649E+04,-3.1911E+03, 7.0180E+03, 1.9862E+03, 5.4146E+01, & 1.4193E+03, 3.5325E+03, 1.2990E+04, 1.9202E+03,-4.1577E+04, &-3.6210E+03, 7.5891E+03, 2.0823E+03, 6.0556E+01, 1.6705E+03, & 3.6979E+03, 1.3664E+04, 1.8183E+03,-4.3643E+04,-4.0505E+03, & 8.1520E+03, 2.1988E+03, 8.8450E+01, 1.8763E+03, 3.9281E+03, & 1.4427E+04, 1.7944E+03,-4.6000E+04,-4.4724E+03, 8.7616E+03, & 2.3413E+03, 1.3066E+02, 2.0515E+03, 4.2101E+03, 1.5254E+04, & 1.8910E+03,-4.8678E+04,-4.8678E+03, 9.4122E+03, 2.5201E+03/ data ((A360(ix,js),js=1,9),ix=21,25)/ & 1.9627E+02, 2.1888E+03, 4.5373E+03, 1.6134E+04, 2.1484E+03, &-5.1718E+04,-5.2303E+03, 1.0080E+04, 2.7544E+03, 2.9253E+02, & 2.2999E+03, 4.9057E+03, 1.7119E+04, 2.5685E+03,-5.5267E+04, &-5.5588E+03, 1.0805E+04, 3.0452E+03, 4.2069E+02, 2.4046E+03, & 5.3221E+03, 1.8287E+04, 3.1724E+03,-5.9480E+04,-5.8650E+03, & 1.1590E+04, 3.4422E+03, 5.6565E+02, 2.5471E+03, 5.7802E+03, & 1.9750E+04, 3.9225E+03,-6.4671E+04,-6.1721E+03, 1.2437E+04, & 4.0778E+03, 6.8164E+02, 2.7280E+03, 6.2513E+03, 2.1474E+04, & 4.7161E+03,-7.0534E+04,-6.4551E+03, 1.3250E+04, 5.0017E+03/ data ((A360(ix,js),js=1,9),ix=26,31)/ & 7.5815E+02, 2.9531E+03, 6.7397E+03, 2.3465E+04, 5.5030E+03, &-7.7007E+04,-6.7321E+03, 1.4003E+04, 6.2391E+03, 8.1943E+02, & 3.2232E+03, 7.2806E+03, 2.5738E+04, 6.2592E+03,-8.4190E+04, &-7.0625E+03, 1.4724E+04, 7.6781E+03, 8.8343E+02, 3.5281E+03, & 7.9013E+03, 2.8287E+04, 6.9149E+03,-9.2081E+04,-7.4700E+03, & 1.5391E+04, 9.1934E+03, 9.6259E+02, 3.8750E+03, 8.6310E+03, & 3.0544E+04, 7.4519E+03,-9.9873E+04,-7.8232E+03, 1.5896E+04, & 1.0647E+04, 1.0657E+03, 4.2669E+03, 9.4391E+03, 3.1950E+04, & 7.2098E+03,-1.0632E+05,-7.8389E+03, 1.6552E+04, 1.1656E+04, & 1.1918E+03, 4.7413E+03, 1.0192E+04, 2.9994E+04, 6.6001E+03, &-1.0795E+05,-7.1220E+03, 1.7021E+04, 1.1715E+04/ data ((A360(ix,js),js=1,9),ix=32,37)/ & 1.3448E+03, 5.3049E+03, 1.0804E+04, 2.8152E+04, 3.0057E+03, &-1.0807E+05,-5.8296E+03, 1.9545E+04, 1.0393E+04, 1.5321E+03, & 6.0282E+03, 1.1182E+04, 2.4855E+04, 3.6060E+03,-1.0871E+05, &-5.2856E+03, 2.0777E+04, 8.2516E+03, 1.7158E+03, 6.8949E+03, & 1.0920E+04, 2.3437E+04, 1.3298E+03,-1.1132E+05,-4.5135E+03, & 2.5203E+04, 6.1970E+03, 1.9016E+03, 7.9667E+03, 9.9618E+03, & 2.0525E+04, 4.2817E+03,-1.1362E+05,-5.0031E+03, 2.7194E+04, & 4.4984E+03, 2.1244E+03, 9.3828E+03, 7.5590E+03, 1.9678E+04, & 7.0220E+03,-1.2033E+05,-4.9126E+03, 3.1851E+04, 3.8717E+03, & 2.3664E+03, 1.1265E+04, 5.6693E+03, 2.1491E+04, 1.1265E+04, &-1.3697E+05,-5.7812E+03, 4.0486E+04, 2.9647E+03/ data ((A360(ix,js),js=1,9),ix=38,43)/ & 2.7913E+03, 1.4030E+04, 4.8102E+03, 2.6050E+04, 1.3722E+04, &-1.7247E+05,-5.1546E+03, 5.9806E+04, 4.0151E+03, 3.2048E+03, & 1.6855E+04, 3.8414E+03, 3.0267E+04, 1.8679E+04,-2.0901E+05, &-5.4977E+03, 7.7718E+04, 5.9062E+03, 3.4514E+03, 1.8244E+04, & 1.8858E+03, 3.4031E+04, 2.2001E+04,-2.3196E+05,-6.1386E+03, & 8.9356E+04, 7.6331E+03, 3.1735E+03, 1.7858E+04, 1.5239E+03, & 3.9874E+04, 3.0144E+04,-2.3702E+05,-1.2064E+04, 8.4326E+04, & 5.9467E+03, 2.8899E+03, 1.7226E+04, 3.6062E+03, 5.0991E+04, & 3.8346E+04,-2.5185E+05,-1.9277E+04, 8.1480E+04, 4.3808E+03, & 2.6308E+03, 1.9510E+04, 6.6321E+03, 5.9251E+04, 3.5217E+04, &-2.6413E+05,-2.5896E+04, 8.2704E+04, 2.2345E+03/ data ((B360(ix,js),js=1,9),ix=1,5)/ & 3.3220E+03, 0.0000E+00, 7.3567E+03, 7.6803E+04, 3.6454E+04, &-2.4080E+05, 4.5492E+03, 1.9812E+04, 2.4305E+04, 3.2356E+03, & 0.0000E+00, 8.9320E+03, 7.6994E+04, 3.5834E+04,-2.4231E+05, & 2.8598E+03, 2.5713E+04, 2.2445E+04, 3.4175E+03, 0.0000E+00, & 1.2232E+04, 7.7599E+04, 2.9999E+04,-2.4245E+05,-8.9565E+02, & 2.5185E+04, 2.2348E+04, 3.8227E+03, 0.0000E+00, 1.8057E+04, & 7.7878E+04, 2.1680E+04,-2.4276E+05,-4.7738E+03, 2.1046E+04, & 2.3454E+04, 4.1639E+03, 0.0000E+00, 2.5061E+04, 7.5590E+04, & 1.6597E+04,-2.4395E+05,-7.3727E+03, 1.8454E+04, 2.4034E+04/ data ((B360(ix,js),js=1,9),ix=6,10)/ & 2.0226E+03, 3.6553E+03, 2.7663E+04, 7.4402E+04, 1.3451E+04, &-2.4500E+05,-9.0677E+03, 1.7600E+04, 2.4039E+04, 2.3371E+03, & 4.9460E+03, 2.6742E+04, 7.3298E+04, 1.1572E+04,-2.4282E+05, &-1.0360E+04, 1.8247E+04, 2.3128E+04, 2.5821E+03, 6.5452E+03, & 2.3738E+04, 7.1130E+04, 9.2331E+03,-2.3376E+05,-1.1453E+04, & 1.8905E+04, 2.1477E+04, 2.8407E+03, 8.7467E+03, 2.1169E+04, & 6.9967E+04, 6.3343E+03,-2.2645E+05,-1.2747E+04, 1.9798E+04, & 2.0046E+04, 3.1777E+03, 1.1663E+04, 2.0303E+04, 7.1744E+04, & 3.4971E+03,-2.2845E+05,-1.4826E+04, 2.1600E+04, 1.9387E+04/ data ((B360(ix,js),js=1,9),ix=11,15)/ & 3.5824E+03, 1.4759E+04, 2.0870E+04, 7.5297E+04, 1.1064E+03, &-2.3739E+05,-1.7845E+04, 2.4720E+04, 1.8972E+04, 4.0208E+03, & 1.7477E+04, 2.2103E+04, 7.8013E+04,-1.3374E+03,-2.4578E+05, &-2.1347E+04, 2.7998E+04, 1.8416E+04, 4.4998E+03, 1.9891E+04, & 2.3925E+04, 8.0263E+04,-4.0759E+03,-2.5531E+05,-2.4750E+04, & 3.1873E+04, 1.7813E+04, 4.9537E+03, 2.1855E+04, 2.5057E+04, & 8.0524E+04,-6.6805E+03,-2.5945E+05,-2.7552E+04, 3.5276E+04, & 1.6945E+04, 5.3966E+03, 2.3717E+04, 2.5315E+04, 7.9668E+04, &-8.9218E+03,-2.5948E+05,-2.9806E+04, 3.8224E+04, 1.5980E+04/ data ((B360(ix,js),js=1,9),ix=16,20)/ & 5.8192E+03, 2.5593E+04, 2.5201E+04, 8.0720E+04,-9.4973E+03, &-2.6599E+05,-3.1836E+04, 4.2603E+04, 1.5458E+04, 6.2188E+03, & 2.7630E+04, 2.4953E+04, 8.2660E+04,-9.0242E+03,-2.7596E+05, &-3.3635E+04, 4.7478E+04, 1.5519E+04, 6.5500E+03, 2.9757E+04, & 2.5013E+04, 8.6318E+04,-7.6888E+03,-2.9146E+05,-3.5697E+04, & 5.3080E+04, 1.6177E+04, 8.8577E+03, 2.6042E+04, 2.7715E+04, & 9.1263E+04,-5.6446E+03,-3.1026E+05,-3.8081E+04, 5.8972E+04, & 1.7414E+04, 1.2119E+04, 2.1474E+04, 3.1224E+04, 9.7933E+04, &-3.2719E+03,-3.3360E+05,-4.1185E+04, 6.4976E+04, 1.9369E+04/ data ((B360(ix,js),js=1,9),ix=21,25)/ & 1.6921E+04, 1.6598E+04, 3.5627E+04, 1.0752E+05,-5.1939E+02, &-3.6518E+05,-4.5432E+04, 7.1580E+04, 2.2536E+04, 2.3498E+04, & 1.2022E+04, 4.0689E+04, 1.2032E+05, 2.3753E+03,-4.0492E+05, &-5.0818E+04, 7.8740E+04, 2.7036E+04, 3.1012E+04, 8.2263E+03, & 4.6120E+04, 1.3578E+05, 5.0895E+03,-4.4899E+05,-5.7140E+04, & 8.5635E+04, 3.2388E+04, 3.6715E+04, 5.6657E+03, 5.0954E+04, & 1.5158E+05, 6.0793E+03,-4.8814E+05,-6.2938E+04, 8.8307E+04, & 3.9590E+04, 4.0345E+04, 4.5077E+03, 5.6889E+04, 1.7253E+05, & 4.4000E+03,-5.3621E+05,-7.0080E+04, 8.9426E+04, 5.0068E+04/ data ((B360(ix,js),js=1,9),ix=26,31)/ & 4.3147E+04, 4.5945E+03, 6.5637E+04, 2.0208E+05,-1.4353E+03, &-6.0197E+05,-8.0122E+04, 9.0303E+04, 6.4638E+04, 4.5803E+04, & 6.1820E+03, 7.8481E+04, 2.4208E+05,-1.2292E+04,-6.9085E+05, &-9.4005E+04, 9.2021E+04, 8.2655E+04, 4.8930E+04, 9.7251E+03, & 9.6788E+04, 2.9438E+05,-2.8855E+04,-8.0949E+05,-1.1269E+05, & 9.5309E+04, 1.0409E+05, 5.2153E+04, 1.5736E+04, 1.2112E+05, & 3.5047E+05,-4.9630E+04,-9.4867E+05,-1.3405E+05, 9.9646E+04, & 1.2758E+05, 5.5289E+04, 2.5000E+04, 1.5248E+05, 4.0443E+05, &-7.8986E+04,-1.1015E+06,-1.5585E+05, 1.0908E+05, 1.4991E+05, & 5.8421E+04, 3.7738E+04, 1.8933E+05, 4.1350E+05,-1.1087E+05, &-1.2104E+06,-1.6912E+05, 1.1441E+05, 1.6453E+05/ data ((B360(ix,js),js=1,9),ix=32,37)/ & 6.2520E+04, 5.6218E+04, 2.3453E+05, 4.2779E+05,-1.7334E+05, &-1.3347E+06,-1.8177E+05, 1.4998E+05, 1.6136E+05, 6.7847E+04, & 8.0624E+04, 2.7957E+05, 3.9240E+05,-1.8747E+05,-1.4391E+06, &-1.9884E+05, 1.7445E+05, 1.3890E+05, 7.4315E+04, 1.1285E+05, & 3.1347E+05, 3.8243E+05,-2.4294E+05,-1.5620E+06,-2.2047E+05, & 2.3239E+05, 1.1589E+05, 8.1983E+04, 1.5529E+05, 3.2509E+05, & 3.1830E+05,-2.2382E+05,-1.6485E+06,-2.4782E+05, 2.6696E+05, & 9.0607E+04, 9.1461E+04, 2.1130E+05, 2.8260E+05, 2.8233E+05, &-2.0994E+05,-1.7623E+06,-2.6761E+05, 3.2415E+05, 7.9589E+04, & 1.0184E+05, 2.8912E+05, 2.4285E+05, 2.8782E+05,-1.9263E+05, &-2.0088E+06,-3.0816E+05, 4.3060E+05, 6.2821E+04/ data ((B360(ix,js),js=1,9),ix=38,43)/ & 1.1806E+05, 4.0223E+05, 2.2517E+05, 3.3973E+05,-2.3078E+05, &-2.5372E+06,-3.5716E+05, 6.5288E+05, 7.7642E+04, 1.3586E+05, & 5.2631E+05, 1.8982E+05, 3.9760E+05,-2.1883E+05,-3.1397E+06, &-4.0467E+05, 8.8836E+05, 1.0255E+05, 1.5133E+05, 6.1562E+05, & 1.0783E+05, 4.8241E+05,-1.6362E+05,-3.6844E+06,-4.1372E+05, & 1.1113E+06, 1.2212E+05, 1.5195E+05, 6.5082E+05, 7.2539E+04, & 6.3688E+05, 5.3880E+04,-4.1600E+06,-4.5347E+05, 1.2093E+06, & 9.2240E+04, 1.4336E+05, 6.1671E+05, 1.2423E+05, 8.1556E+05, & 1.3919E+05,-4.2989E+06,-5.6767E+05, 1.1435E+06, 6.9843E+04, & 1.3071E+05, 6.5514E+05, 1.9799E+05, 8.8504E+05, 9.0991E+04, &-4.2280E+06,-6.1482E+05, 1.0996E+06, 3.0361E+04/ data ((A540(ix,js),js=1,9),ix=1,5)/ & 1.4046E+01, 0.0000E+00, 8.4121E+01, 2.7000E+03, 4.9855E+03, &-1.4428E+04, 1.1932E+03, 3.1112E+03, 5.2292E+02, 1.8382E+01, & 0.0000E+00, 1.3808E+02, 3.5508E+03, 5.2498E+03,-1.7036E+04, & 1.1096E+03, 3.6461E+03, 7.1571E+02, 2.2854E+01, 0.0000E+00, & 2.2670E+02, 4.5327E+03, 5.5670E+03,-2.0072E+04, 9.0698E+02, & 4.0662E+03, 9.6773E+02, 2.6788E+01, 0.0000E+00, 3.6639E+02, & 5.5619E+03, 5.8004E+03,-2.3126E+04, 6.7679E+02, 4.2755E+03, & 1.2365E+03, 3.1310E+01, 0.0000E+00, 5.8732E+02, 6.6192E+03, & 5.8620E+03,-2.5990E+04, 3.3715E+02, 4.3206E+03, 1.5084E+03/ data ((A540(ix,js),js=1,9),ix=6,10)/ & 6.4043E+00, 4.6835E+01, 8.5364E+02, 7.7821E+03, 5.6476E+03, &-2.8674E+04,-3.2399E+01, 4.3542E+03, 1.7554E+03, 9.6264E+00, & 7.8879E+01, 1.1548E+03, 8.8791E+03, 5.2330E+03,-3.1025E+04, &-3.5404E+02, 4.4739E+03, 1.9480E+03, 1.3550E+01, 1.3075E+02, & 1.4391E+03, 9.8097E+03, 4.8032E+03,-3.2995E+04,-6.8661E+02, & 4.6960E+03, 2.0777E+03, 1.7714E+01, 2.0669E+02, 1.6861E+03, & 1.0641E+04, 4.3793E+03,-3.4751E+04,-1.0292E+03, 5.0045E+03, & 2.1496E+03, 2.1444E+01, 3.0070E+02, 1.9210E+03, 1.1366E+04, & 3.9995E+03,-3.6335E+04,-1.3571E+03, 5.3448E+03, 2.1701E+03/ data ((A540(ix,js),js=1,9),ix=11,15)/ & 2.4860E+01, 4.0590E+02, 2.1813E+03, 1.1988E+04, 3.7302E+03, &-3.7857E+04,-1.6700E+03, 5.7183E+03, 2.1533E+03, 2.8308E+01, & 5.2176E+02, 2.5305E+03, 1.2550E+04, 3.5359E+03,-3.9488E+04, &-2.0723E+03, 6.1238E+03, 2.1332E+03, 3.2231E+01, 6.6569E+02, & 2.9583E+03, 1.3080E+04, 3.2613E+03,-4.1264E+04,-2.4406E+03, & 6.6123E+03, 2.1193E+03, 3.6139E+01, 8.3147E+02, 3.3614E+03, & 1.3578E+04, 3.0508E+03,-4.3167E+04,-2.8085E+03, 7.2030E+03, & 2.1301E+03, 4.0256E+01, 1.0297E+03, 3.7022E+03, 1.4171E+04, & 2.9128E+03,-4.5336E+04,-3.2104E+03, 7.8946E+03, 2.1776E+03/ data ((A540(ix,js),js=1,9),ix=16,20)/ & 4.4889E+01, 1.2599E+03, 3.9615E+03, 1.4830E+04, 2.8030E+03, &-4.7625E+04,-3.6641E+03, 8.6225E+03, 2.2558E+03, 5.0015E+01, & 1.5224E+03, 4.1757E+03, 1.5557E+04, 2.7204E+03,-5.0044E+04, &-4.1690E+03, 9.3569E+03, 2.3618E+03, 5.6369E+01, 1.8214E+03, & 4.3850E+03, 1.6354E+04, 2.6560E+03,-5.2643E+04,-4.6917E+03, & 1.0073E+04, 2.4979E+03, 8.3103E+01, 2.0790E+03, 4.6427E+03, & 1.7254E+04, 2.6678E+03,-5.5501E+04,-5.2101E+03, 1.0809E+04, & 2.6617E+03, 1.2393E+02, 2.3045E+03, 4.9515E+03, 1.8262E+04, & 2.8040E+03,-5.8738E+04,-5.7019E+03, 1.1579E+04, 2.8680E+03/ data ((A540(ix,js),js=1,9),ix=21,25)/ & 1.8908E+02, 2.4870E+03, 5.3135E+03, 1.9357E+04, 3.1110E+03, &-6.2419E+04,-6.1562E+03, 1.2371E+04, 3.1342E+03, 2.8397E+02, & 2.6277E+03, 5.7304E+03, 2.0564E+04, 3.6012E+03,-6.6647E+04, &-6.5563E+03, 1.3212E+04, 3.4531E+03, 4.1265E+02, 2.7604E+03, & 6.2057E+03, 2.1959E+04, 4.3164E+03,-7.1623E+04,-6.9163E+03, & 1.4129E+04, 3.8711E+03, 5.5908E+02, 2.9304E+03, 6.7324E+03, & 2.3661E+04, 5.2276E+03,-7.7697E+04,-7.2578E+03, 1.5148E+04, & 4.5144E+03, 6.7735E+02, 3.1384E+03, 7.2790E+03, 2.5647E+04, & 6.2256E+03,-8.4578E+04,-7.5573E+03, 1.6206E+04, 5.4499E+03/ data ((A540(ix,js),js=1,9),ix=26,31)/ & 7.6095E+02, 3.3936E+03, 7.8572E+03, 2.7951E+04, 7.2729E+03, &-9.2312E+04,-7.8460E+03, 1.7305E+04, 6.7468E+03, 8.2152E+02, & 3.6897E+03, 8.4910E+03, 3.0565E+04, 8.3712E+03,-1.0096E+05, &-8.1788E+03, 1.8426E+04, 8.3335E+03, 8.9734E+02, 4.0371E+03, & 9.2179E+03, 3.3541E+04, 9.4243E+03,-1.1068E+05,-8.6295E+03, & 1.9539E+04, 1.0194E+04, 9.8654E+02, 4.4310E+03, 1.0051E+04, & 3.6202E+04, 1.0341E+04,-1.2033E+05,-9.0465E+03, 2.0323E+04, & 1.2170E+04, 1.1010E+03, 4.8898E+03, 1.0979E+04, 3.8051E+04, & 1.0183E+04,-1.2865E+05,-9.1689E+03, 2.1141E+04, 1.3739E+04, & 1.2523E+03, 5.4507E+03, 1.1845E+04, 3.6056E+04, 9.2641E+03, &-1.3103E+05,-8.5702E+03, 2.1465E+04, 1.4032E+04/ data ((A540(ix,js),js=1,9),ix=32,37)/ & 1.4290E+03, 6.1493E+03, 1.2630E+04, 3.4660E+04, 4.6113E+03, &-1.3253E+05,-7.6923E+03, 2.5385E+04, 1.1814E+04, 1.6393E+03, & 6.9603E+03, 1.3045E+04, 3.1292E+04, 5.0250E+03,-1.3241E+05, &-8.0218E+03, 2.7248E+04, 8.1799E+03, 1.8422E+03, 7.9125E+03, & 1.2850E+04, 2.9919E+04, 1.6539E+03,-1.3374E+05,-7.5523E+03, & 3.2067E+04, 5.4650E+03, 2.0556E+03, 9.1198E+03, 1.1952E+04, & 2.6229E+04, 5.4064E+03,-1.3513E+05,-8.5919E+03, 3.3805E+04, & 3.3909E+03, 2.3195E+03, 1.0765E+04, 9.4232E+03, 2.4680E+04, & 7.8282E+03,-1.4260E+05,-8.1716E+03, 3.9314E+04, 3.2746E+03, & 2.6741E+03, 1.3056E+04, 7.3683E+03, 2.5901E+04, 1.1853E+04, &-1.6172E+05,-8.6679E+03, 4.9358E+04, 3.1283E+03/ data ((A540(ix,js),js=1,9),ix=38,43)/ & 3.1274E+03, 1.6190E+04, 6.6171E+03, 2.9852E+04, 1.3996E+04, &-2.0033E+05,-7.1656E+03, 7.0586E+04, 5.0072E+03, 3.5342E+03, & 1.9357E+04, 5.6559E+03, 3.3611E+04, 1.7828E+04,-2.3702E+05, &-6.2284E+03, 8.8163E+04, 7.3170E+03, 3.7834E+03, 2.0965E+04, & 3.2588E+03, 3.7010E+04, 2.1189E+04,-2.5903E+05,-5.8621E+03, & 9.8554E+04, 8.6223E+03, 3.5046E+03, 2.0797E+04, 1.9896E+03, & 4.3329E+04, 3.0815E+04,-2.6553E+05,-1.1030E+04, 9.3433E+04, & 6.0978E+03, 3.3389E+03, 2.0607E+04, 2.5567E+03, 5.6928E+04, & 4.3500E+04,-2.8735E+05,-1.7834E+04, 8.9701E+04, 3.7020E+03, & 3.2675E+03, 2.3717E+04, 3.5488E+03, 7.1174E+04, 4.9418E+04, &-3.1866E+05,-2.5099E+04, 9.4302E+04, 1.5043E+03/ data ((B540(ix,js),js=1,9),ix=1,5)/ & 2.9709E+03, 0.0000E+00, 7.2631E+03, 8.2400E+04, 4.6777E+04, &-2.6896E+05, 8.8736E+03, 2.2636E+04, 2.6978E+04, 3.1238E+03, & 0.0000E+00, 9.6559E+03, 8.7263E+04, 4.5581E+04,-2.8113E+05, & 6.5338E+03, 3.0570E+04, 2.5612E+04, 3.4594E+03, 0.0000E+00, & 1.3781E+04, 9.1771E+04, 4.2054E+04,-2.9581E+05, 2.4720E+03, & 3.2871E+04, 2.6461E+04, 3.9369E+03, 0.0000E+00, 2.0602E+04, & 9.1745E+04, 3.0580E+04,-2.9239E+05,-2.5695E+03, 2.6492E+04, & 2.8006E+04, 4.3495E+03, 0.0000E+00, 2.8814E+04, 8.7356E+04, & 2.1482E+04,-2.8467E+05,-6.3525E+03, 2.1028E+04, 2.8623E+04/ data ((B540(ix,js),js=1,9),ix=6,10)/ & 2.0005E+03, 4.0531E+03, 3.1958E+04, 8.4074E+04, 1.6201E+04, &-2.7909E+05,-8.6225E+03, 1.8720E+04, 2.8409E+04, 2.3629E+03, & 5.5993E+03, 3.1028E+04, 8.2417E+04, 1.3922E+04,-2.7599E+05, &-1.0282E+04, 1.9284E+04, 2.7425E+04, 2.6727E+03, 7.5414E+03, & 2.7733E+04, 8.1705E+04, 1.1962E+04,-2.7148E+05,-1.1901E+04, & 2.0897E+04, 2.5922E+04, 2.9924E+03, 1.0195E+04, 2.4967E+04, & 8.2254E+04, 9.0794E+03,-2.6878E+05,-1.3758E+04, 2.2862E+04, & 2.4609E+04, 3.3628E+03, 1.3635E+04, 2.4057E+04, 8.5280E+04, & 5.7370E+03,-2.7339E+05,-1.6337E+04, 2.5177E+04, 2.4060E+04/ data ((B540(ix,js),js=1,9),ix=11,15)/ & 3.7910E+03, 1.7353E+04, 2.4740E+04, 8.9874E+04, 2.8386E+03, &-2.8470E+05,-2.0176E+04, 2.9246E+04, 2.3552E+04, 4.2537E+03, & 2.0528E+04, 2.6524E+04, 9.4192E+04,-1.9538E+02,-2.9713E+05, &-2.4840E+04, 3.3537E+04, 2.3038E+04, 4.7553E+03, 2.3449E+04, & 2.8849E+04, 9.7935E+04,-3.3681E+03,-3.1112E+05,-2.9373E+04, & 3.8841E+04, 2.2236E+04, 5.2435E+03, 2.5812E+04, 3.0777E+04, & 1.0004E+05,-6.4184E+03,-3.2110E+05,-3.3387E+04, 4.3876E+04, & 2.1211E+04, 5.7016E+03, 2.8053E+04, 3.1613E+04, 1.0025E+05, &-9.2545E+03,-3.2451E+05,-3.6675E+04, 4.8014E+04, 1.9899E+04/ data ((B540(ix,js),js=1,9),ix=16,20)/ & 6.1442E+03, 3.0360E+04, 3.1866E+04, 1.0151E+05,-1.0762E+04, &-3.3154E+05,-3.9524E+04, 5.3002E+04, 1.8817E+04, 6.5885E+03, & 3.3054E+04, 3.1916E+04, 1.0396E+05,-1.0604E+04,-3.4429E+05, &-4.1917E+04, 5.8965E+04, 1.8358E+04, 6.9707E+03, 3.5930E+04, & 3.2193E+04, 1.0817E+05,-8.9122E+03,-3.6374E+05,-4.4279E+04, & 6.5938E+04, 1.8630E+04, 9.5116E+03, 3.2665E+04, 3.5521E+04, & 1.1452E+05,-5.9078E+03,-3.9000E+05,-4.7090E+04, 7.3913E+04, & 1.9693E+04, 1.3035E+04, 2.8272E+04, 3.9666E+04, 1.2272E+05, &-1.9538E+03,-4.2149E+05,-5.0485E+04, 8.2204E+04, 2.1617E+04/ data ((B540(ix,js),js=1,9),ix=21,25)/ & 1.8173E+04, 2.3394E+04, 4.4584E+04, 1.3380E+05, 2.9300E+03, &-4.6142E+05,-5.4863E+04, 9.1181E+04, 2.4808E+04, 2.5309E+04, & 1.8719E+04, 5.0230E+04, 1.4885E+05, 8.6209E+03,-5.1236E+05, &-6.0511E+04, 1.0139E+05, 2.9436E+04, 3.3629E+04, 1.4750E+04, & 5.6231E+04, 1.6736E+05, 1.4927E+04,-5.7045E+05,-6.7277E+04, & 1.1217E+05, 3.5060E+04, 3.9981E+04, 1.1960E+04, 6.1755E+04, & 1.8744E+05, 2.0137E+04,-6.2674E+05,-7.3834E+04, 1.1934E+05, & 4.3078E+04, 4.4113E+04, 1.0697E+04, 6.8456E+04, 2.1400E+05, & 2.3263E+04,-6.9595E+05,-8.2136E+04, 1.2566E+05, 5.5471E+04/ data ((B540(ix,js),js=1,9),ix=26,31)/ & 4.7268E+04, 1.0742E+04, 7.7893E+04, 2.5024E+05, 2.2431E+04, &-7.8521E+05,-9.3714E+04, 1.3164E+05, 7.3653E+04, 5.0261E+04, & 1.2320E+04, 9.1514E+04, 2.9871E+05, 1.5721E+04,-9.0003E+05, &-1.0984E+05, 1.3732E+05, 9.7259E+04, 5.3382E+04, 1.5997E+04, & 1.1025E+05, 3.6030E+05, 1.6491E+03,-1.0428E+06,-1.3104E+05, & 1.4263E+05, 1.2572E+05, 5.6536E+04, 2.2075E+04, 1.3502E+05, & 4.2540E+05,-1.9356E+04,-1.2017E+06,-1.5497E+05, 1.4620E+05, & 1.5667E+05, 5.9668E+04, 3.1332E+04, 1.6712E+05, 4.8760E+05, &-5.2497E+04,-1.3721E+06,-1.7921E+05, 1.5491E+05, 1.8500E+05, & 6.3128E+04, 4.4306E+04, 2.0538E+05, 5.0036E+05,-8.8315E+04, &-1.5009E+06,-1.9544E+05, 1.6286E+05, 2.0060E+05/ data ((B540(ix,js),js=1,9),ix=32,37)/ & 6.7252E+04, 6.1350E+04, 2.4981E+05, 5.1866E+05,-1.7393E+05, &-1.6079E+06,-2.1616E+05, 2.0514E+05, 1.8421E+05, 7.2523E+04, & 8.7077E+04, 3.0032E+05, 4.9795E+05,-1.9374E+05,-1.7488E+06, &-2.5154E+05, 2.5270E+05, 1.4108E+05, 7.8557E+04, 1.2085E+05, & 3.4369E+05, 5.0559E+05,-2.7072E+05,-1.8995E+06,-2.8633E+05, & 3.2697E+05, 1.0954E+05, 8.6092E+04, 1.6790E+05, 3.6976E+05, & 4.3787E+05,-2.5173E+05,-2.0235E+06,-3.3252E+05, 3.7066E+05, & 7.9724E+04, 9.6434E+04, 2.3315E+05, 3.3902E+05, 3.9461E+05, &-2.5455E+05,-2.2059E+06,-3.6308E+05, 4.4700E+05, 8.0952E+04, & 1.0946E+05, 3.2289E+05, 3.0583E+05, 3.8625E+05,-2.5579E+05, &-2.5156E+06,-4.1451E+05, 5.7253E+05, 7.8418E+04/ data ((B540(ix,js),js=1,9),ix=38,43)/ & 1.2692E+05, 4.4886E+05, 3.0244E+05, 4.2615E+05,-3.1844E+05, &-3.1255E+06,-4.6924E+05, 8.2340E+05, 1.1046E+05, 1.4441E+05, & 5.8243E+05, 2.7541E+05, 4.7182E+05,-3.3649E+05,-3.7185E+06, &-5.1067E+05, 1.0440E+06, 1.4277E+05, 1.5916E+05, 6.7402E+05, & 1.8070E+05, 5.5047E+05,-2.9324E+05,-4.2164E+06,-5.1434E+05, & 1.2427E+06, 1.5864E+05, 1.6012E+05, 7.2380E+05, 1.1475E+05, & 7.1511E+05,-6.9798E+04,-4.6914E+06,-5.4151E+05, 1.3302E+06, & 1.1012E+05, 1.5471E+05, 7.0419E+05, 1.1310E+05, 8.8568E+05, & 2.4896E+04,-4.6768E+06,-6.3312E+05, 1.1713E+06, 6.4977E+04, & 1.5021E+05, 7.9025E+05, 1.1779E+05, 1.0507E+06, 1.1771E+05, &-4.9507E+06,-6.8570E+05, 1.1955E+06, 2.3632E+04/ data ((A720(ix,js),js=1,9),ix=1,5)/ & 1.1000E+01, 0.0000E+00, 8.1568E+01, 2.8044E+03, 5.3898E+03, &-1.5286E+04, 1.3891E+03, 3.3052E+03, 5.0976E+02, 1.5350E+01, & 0.0000E+00, 1.3680E+02, 3.6955E+03, 5.7560E+03,-1.8147E+04, & 1.3596E+03, 3.9557E+03, 6.9200E+02, 1.9863E+01, 0.0000E+00, & 2.2561E+02, 4.7049E+03, 6.2084E+03,-2.1486E+04, 1.1961E+03, & 4.5137E+03, 9.4594E+02, 2.4344E+01, 0.0000E+00, 3.6564E+02, & 5.7723E+03, 6.6036E+03,-2.4949E+04, 9.9205E+02, 4.8226E+03, & 1.2417E+03, 2.8832E+01, 0.0000E+00, 5.8203E+02, 6.9163E+03, & 6.8325E+03,-2.8322E+04, 6.4255E+02, 4.9144E+03, 1.5622E+03/ data ((A720(ix,js),js=1,9),ix=6,10)/ & 4.8883E+00, 4.4839E+01, 8.4544E+02, 8.2527E+03, 6.7546E+03, &-3.1641E+04, 2.3650E+02, 4.9722E+03, 1.8764E+03, 7.5913E+00, & 7.6632E+01, 1.1494E+03, 9.5596E+03, 6.3700E+03,-3.4538E+04, &-1.2944E+02, 5.1143E+03, 2.1263E+03, 1.1092E+01, 1.2916E+02, & 1.4465E+03, 1.0700E+04, 5.9160E+03,-3.6968E+04,-5.1496E+02, & 5.3844E+03, 2.2995E+03, 1.5051E+01, 2.0677E+02, 1.7088E+03, & 1.1710E+04, 5.4364E+03,-3.9056E+04,-9.2065E+02, 5.7632E+03, & 2.3914E+03, 1.8844E+01, 3.0261E+02, 1.9641E+03, 1.2606E+04, & 4.9993E+03,-4.0960E+04,-1.3152E+03, 6.1798E+03, 2.4231E+03/ data ((A720(ix,js),js=1,9),ix=11,15)/ & 2.2391E+01, 4.0801E+02, 2.2543E+03, 1.3383E+04, 4.6738E+03, &-4.2763E+04,-1.6992E+03, 6.6211E+03, 2.4071E+03, 2.5971E+01, & 5.2138E+02, 2.6560E+03, 1.4114E+04, 4.4287E+03,-4.4726E+04, &-2.1995E+03, 7.0889E+03, 2.3855E+03, 2.9763E+01, 6.6460E+02, & 3.1664E+03, 1.4809E+04, 4.0645E+03,-4.6843E+04,-2.6537E+03, & 7.6437E+03, 2.3613E+03, 3.3666E+01, 8.3830E+02, 3.6616E+03, & 1.5430E+04, 3.7727E+03,-4.9054E+04,-3.0909E+03, 8.3082E+03, & 2.3574E+03, 3.7748E+01, 1.0512E+03, 4.0875E+03, 1.6120E+04, & 3.5818E+03,-5.1529E+04,-3.5443E+03, 9.0933E+03, 2.3920E+03/ data ((A720(ix,js),js=1,9),ix=16,20)/ & 4.2296E+01, 1.3036E+03, 4.4196E+03, 1.6881E+04, 3.4671E+03, &-5.4229E+04,-4.0449E+03, 9.9546E+03, 2.4679E+03, 4.7480E+01, & 1.5965E+03, 4.6912E+03, 1.7692E+04, 3.4059E+03,-5.7053E+04, &-4.6035E+03, 1.0823E+04, 2.5794E+03, 5.3627E+01, 1.9302E+03, & 4.9402E+03, 1.8564E+04, 3.3780E+03,-6.0042E+04,-5.1911E+03, & 1.1665E+04, 2.7244E+03, 7.9167E+01, 2.2283E+03, 5.2269E+03, & 1.9564E+04, 3.4343E+03,-6.3324E+04,-5.7858E+03, 1.2522E+04, & 2.9022E+03, 1.1948E+02, 2.4959E+03, 5.5584E+03, 2.0703E+04, & 3.6185E+03,-6.7027E+04,-6.3585E+03, 1.3402E+04, 3.1300E+03/ data ((A720(ix,js),js=1,9),ix=21,25)/ & 1.8391E+02, 2.7146E+03, 5.9441E+03, 2.1967E+04, 3.9800E+03, &-7.1233E+04,-6.8930E+03, 1.4300E+04, 3.4242E+03, 2.7904E+02, & 2.8863E+03, 6.3937E+03, 2.3372E+04, 4.5327E+03,-7.6066E+04, &-7.3642E+03, 1.5250E+04, 3.7769E+03, 4.0911E+02, 3.0443E+03, & 6.9146E+03, 2.4981E+04, 5.3338E+03,-8.1711E+04,-7.7825E+03, & 1.6272E+04, 4.2287E+03, 5.5647E+02, 3.2369E+03, 7.4969E+03, & 2.6902E+04, 6.3584E+03,-8.8522E+04,-8.1632E+03, 1.7400E+04, & 4.8989E+03, 6.7746E+02, 3.4639E+03, 8.1044E+03, 2.9116E+04, & 7.4937E+03,-9.6201E+04,-8.4797E+03, 1.8601E+04, 5.8474E+03/ data ((A720(ix,js),js=1,9),ix=26,31)/ & 7.6419E+02, 3.7351E+03, 8.7528E+03, 3.1669E+04, 8.7109E+03, &-1.0484E+05,-8.7688E+03, 1.9909E+04, 7.1538E+03, 8.3217E+02, & 4.0555E+03, 9.4681E+03, 3.4560E+04, 1.0047E+04,-1.1462E+05, &-9.0931E+03, 2.1351E+04, 8.7935E+03, 9.1301E+02, 4.4297E+03, & 1.0290E+04, 3.7864E+04, 1.1445E+04,-1.2584E+05,-9.5391E+03, & 2.2897E+04, 1.0786E+04, 1.0164E+03, 4.8721E+03, 1.1249E+04, & 4.0878E+04, 1.2767E+04,-1.3731E+05,-9.9925E+03, 2.4149E+04, & 1.3084E+04, 1.1415E+03, 5.3815E+03, 1.2293E+04, 4.2941E+04, & 1.2877E+04,-1.4705E+05,-1.0110E+04, 2.5288E+04, 1.5084E+04, & 1.2975E+03, 5.9819E+03, 1.3226E+04, 4.0585E+04, 1.1739E+04, &-1.4918E+05,-9.3569E+03, 2.5526E+04, 1.5795E+04/ data ((A720(ix,js),js=1,9),ix=32,37)/ & 1.4828E+03, 6.7270E+03, 1.4083E+04, 3.9252E+04, 6.2124E+03, &-1.5042E+05,-8.4580E+03, 2.9972E+04, 1.3334E+04, 1.7053E+03, & 7.6200E+03, 1.4682E+04, 3.6183E+04, 5.9149E+03,-1.5116E+05, &-9.3969E+03, 3.2817E+04, 9.0065E+03, 1.9276E+03, 8.6904E+03, & 1.4523E+04, 3.5751E+04, 2.4091E+03,-1.5492E+05,-1.0156E+04, & 3.9419E+04, 4.9338E+03, 2.1451E+03, 1.0040E+04, 1.3698E+04, & 3.1820E+04, 7.1113E+03,-1.5760E+05,-1.1397E+04, 3.9737E+04, & 3.4762E+03, 2.4469E+03, 1.1859E+04, 1.0944E+04, 2.9494E+04, & 8.8440E+03,-1.6440E+05,-9.9882E+03, 4.4020E+04, 4.6728E+03, & 2.8411E+03, 1.4376E+04, 8.8420E+03, 3.0261E+04, 1.2589E+04, &-1.8408E+05,-9.3312E+03, 5.2682E+04, 6.4506E+03/ data ((A720(ix,js),js=1,9),ix=38,43)/ & 3.3274E+03, 1.7805E+04, 8.3315E+03, 3.4132E+04, 1.3953E+04, &-2.2590E+05,-7.2006E+03, 7.6708E+04, 8.4541E+03, 3.7593E+03, & 2.1296E+04, 7.4382E+03, 3.7788E+04, 1.7525E+04,-2.6251E+05, &-6.7751E+03, 9.5562E+04, 9.7409E+03, 4.0215E+03, 2.3146E+04, & 4.4673E+03, 3.9980E+04, 1.9665E+04,-2.8011E+05,-6.2491E+03, & 1.0617E+05, 9.7124E+03, 3.7615E+03, 2.3278E+04, 2.4201E+03, & 4.5810E+04, 2.9998E+04,-2.8413E+05,-1.1164E+04, 9.8953E+04, & 6.3735E+03, 3.7240E+03, 2.3588E+04, 2.0106E+03, 6.0596E+04, & 4.6732E+04,-3.1152E+05,-1.7614E+04, 9.4606E+04, 3.4489E+03, & 3.7961E+03, 2.6755E+04, 2.3096E+03, 7.7763E+04, 5.7709E+04, &-3.5374E+05,-2.4401E+04, 1.0057E+05, 8.7498E+02/ data ((B720(ix,js),js=1,9),ix=1,5)/ & 2.8378E+03, 0.0000E+00, 7.4867E+03, 8.7260E+04, 4.9146E+04, &-2.8306E+05, 1.0229E+04, 2.1829E+04, 2.9344E+04, 3.0405E+03, & 0.0000E+00, 1.0134E+04, 9.3386E+04, 4.9056E+04,-2.9980E+05, & 8.3963E+03, 3.1439E+04, 2.7751E+04, 3.4447E+03, 0.0000E+00, & 1.4768E+04, 9.9677E+04, 4.7199E+04,-3.2228E+05, 4.6847E+03, & 3.5809E+04, 2.8838E+04, 3.9763E+03, 0.0000E+00, 2.2241E+04, & 1.0162E+05, 3.8195E+04,-3.2951E+05,-3.5761E+02, 3.1367E+04, & 3.1066E+04, 4.4566E+03, 0.0000E+00, 3.1280E+04, 9.7135E+04, & 2.7683E+04,-3.2113E+05,-4.8036E+03, 2.4666E+04, 3.2099E+04/ data ((B720(ix,js),js=1,9),ix=6,10)/ & 1.9920E+03, 4.3042E+03, 3.4960E+04, 9.2325E+04, 1.9550E+04, &-3.0863E+05,-7.8736E+03, 2.0414E+04, 3.1867E+04, 2.3607E+03, & 6.0047E+03, 3.4158E+04, 9.0097E+04, 1.6570E+04,-3.0384E+05, &-9.9396E+03, 2.0701E+04, 3.0825E+04, 2.7001E+03, 8.1576E+03, & 3.0612E+04, 8.9617E+04, 1.4401E+04,-2.9982E+05,-1.1977E+04, & 2.2595E+04, 2.9273E+04, 3.0360E+03, 1.1081E+04, 2.7697E+04, & 9.1658E+04, 1.1561E+04,-3.0089E+05,-1.4403E+04, 2.5524E+04, & 2.8020E+04, 3.4215E+03, 1.4885E+04, 2.6879E+04, 9.6042E+04, & 7.9992E+03,-3.0885E+05,-1.7495E+04, 2.8632E+04, 2.7565E+04/ data ((B720(ix,js),js=1,9),ix=11,15)/ & 3.8631E+03, 1.9005E+04, 2.7888E+04, 1.0225E+05, 4.8521E+03, &-3.2437E+05,-2.2097E+04, 3.3740E+04, 2.7150E+04, 4.3522E+03, & 2.2605E+04, 3.0259E+04, 1.0826E+05, 1.4579E+03,-3.4145E+05, &-2.7813E+04, 3.9189E+04, 2.6718E+04, 4.8734E+03, 2.5856E+04, & 3.3158E+04, 1.1340E+05,-2.1817E+03,-3.5946E+05,-3.3255E+04, & 4.5530E+04, 2.5944E+04, 5.3812E+03, 2.8555E+04, 3.5608E+04, & 1.1671E+05,-5.6934E+03,-3.7309E+05,-3.8155E+04, 5.1594E+04, & 2.4865E+04, 5.8581E+03, 3.1107E+04, 3.6904E+04, 1.1784E+05, &-9.1119E+03,-3.7911E+05,-4.2332E+04, 5.6614E+04, 2.3383E+04/ data ((B720(ix,js),js=1,9),ix=16,20)/ & 6.3262E+03, 3.3720E+04, 3.7274E+04, 1.1891E+05,-1.1690E+04, &-3.8488E+05,-4.5891E+04, 6.1655E+04, 2.1870E+04, 6.8035E+03, & 3.6904E+04, 3.7600E+04, 1.2174E+05,-1.2075E+04,-3.9908E+05, &-4.8826E+04, 6.8095E+04, 2.1052E+04, 7.2152E+03, 4.0328E+04, & 3.8374E+04, 1.2712E+05,-1.0373E+04,-4.2320E+05,-5.1740E+04, & 7.6305E+04, 2.0951E+04, 9.8640E+03, 3.7426E+04, 4.2372E+04, & 1.3471E+05,-6.8135E+03,-4.5532E+05,-5.4983E+04, 8.5856E+04, & 2.1689E+04, 1.3560E+04, 3.3322E+04, 4.7176E+04, 1.4419E+05, &-1.7431E+03,-4.9373E+05,-5.8662E+04, 9.5890E+04, 2.3540E+04/ data ((B720(ix,js),js=1,9),ix=21,25)/ & 1.8962E+04, 2.8655E+04, 5.2627E+04, 1.5648E+05, 4.6977E+03, &-5.4094E+05,-6.3178E+04, 1.0669E+05, 2.6814E+04, 2.6434E+04, & 2.4076E+04, 5.8634E+04, 1.7285E+05, 1.2578E+04,-5.9982E+05, &-6.8826E+04, 1.1902E+05, 3.1514E+04, 3.5147E+04, 2.0130E+04, & 6.5020E+04, 1.9328E+05, 2.1789E+04,-6.6785E+05,-7.5631E+04, & 1.3244E+05, 3.7235E+04, 4.1827E+04, 1.7320E+04, 7.0904E+04, & 2.1618E+05, 3.0623E+04,-7.3675E+05,-8.2378E+04, 1.4290E+05, & 4.5585E+04, 4.6294E+04, 1.6071E+04, 7.8131E+04, 2.4682E+05, & 3.8078E+04,-8.2274E+05,-9.1157E+04, 1.5361E+05, 5.8726E+04/ data ((B720(ix,js),js=1,9),ix=26,31)/ & 4.9878E+04, 1.6065E+04, 8.8228E+04, 2.8824E+05, 4.2519E+04, &-9.3298E+05,-1.0344E+05, 1.6515E+05, 7.8533E+04, 5.3229E+04, & 1.7740E+04, 1.0255E+05, 3.4331E+05, 4.1675E+04,-1.0732E+06, &-1.2083E+05, 1.7696E+05, 1.0553E+05, 5.6585E+04, 2.1720E+04, & 1.2220E+05, 4.1397E+05, 3.3029E+04,-1.2465E+06,-1.4427E+05, & 1.8802E+05, 1.3975E+05, 5.9316E+04, 2.8288E+04, 1.4733E+05, & 4.8733E+05, 1.5630E+04,-1.4301E+06,-1.7016E+05, 1.9433E+05, & 1.7788E+05, 6.2252E+04, 3.8013E+04, 1.8011E+05, 5.5712E+05, &-1.8402E+04,-1.6194E+06,-1.9622E+05, 2.0377E+05, 2.1356E+05, & 6.5578E+04, 5.1423E+04, 2.1882E+05, 5.6668E+05,-5.8779E+04, &-1.7453E+06,-2.1194E+05, 2.1181E+05, 2.3210E+05/ data ((B720(ix,js),js=1,9),ix=32,37)/ & 6.9606E+04, 6.8196E+04, 2.6312E+05, 5.9063E+05,-1.5744E+05, &-1.8524E+06,-2.3319E+05, 2.6085E+05, 2.1170E+05, 7.5218E+04, & 9.4375E+04, 3.1888E+05, 5.8441E+05,-1.8982E+05,-2.0276E+06, &-2.8095E+05, 3.2730E+05, 1.5852E+05, 8.1818E+04, 1.2971E+05, & 3.6827E+05, 6.1914E+05,-2.8165E+05,-2.2420E+06,-3.4355E+05, & 4.3737E+05, 1.0553E+05, 8.9628E+04, 1.7800E+05, 4.0374E+05, & 5.5218E+05,-2.6089E+05,-2.4102E+06,-4.0173E+05, 4.6432E+05, & 8.6117E+04, 1.0036E+05, 2.4647E+05, 3.7601E+05, 4.9780E+05, &-2.7995E+05,-2.6053E+06,-4.2444E+05, 5.1942E+05, 1.1274E+05, & 1.1449E+05, 3.4157E+05, 3.5045E+05, 4.8110E+05,-2.9572E+05, &-2.9462E+06,-4.6715E+05, 6.1840E+05, 1.4715E+05/ data ((B720(ix,js),js=1,9),ix=38,43)/ & 1.3273E+05, 4.7433E+05, 3.6299E+05, 5.2029E+05,-3.8660E+05, &-3.6364E+06,-5.3004E+05, 9.1163E+05, 1.8527E+05, 1.5094E+05, & 6.1556E+05, 3.4710E+05, 5.6221E+05,-4.2073E+05,-4.2353E+06, &-5.8904E+05, 1.1559E+06, 2.0012E+05, 1.6568E+05, 7.1379E+05, & 2.4118E+05, 6.1735E+05,-3.9948E+05,-4.6347E+06,-5.9224E+05, & 1.3508E+06, 1.8977E+05, 1.6599E+05, 7.6711E+05, 1.5253E+05, & 7.5990E+05,-2.0122E+05,-4.9316E+06,-6.2692E+05, 1.3721E+06, & 1.2698E+05, 1.6400E+05, 7.6327E+05, 1.1905E+05, 9.5620E+05, &-5.0089E+04,-5.0316E+06,-7.0061E+05, 1.2124E+06, 6.6342E+04, & 1.6539E+05, 8.6958E+05, 9.0836E+04, 1.1753E+06, 8.8077E+04, &-5.4751E+06,-7.6386E+05, 1.2589E+06, 2.1098E+04/ end block data