diff --git a/plasim/src/radmod.f90 b/plasim/src/radmod.f90 index cf62925..e971c7e 100644 --- a/plasim/src/radmod.f90 +++ b/plasim/src/radmod.f90 @@ -16,11 +16,15 @@ ! parameter(SBK = 5.67E-8) ! Stefan-Bolzman Const. +! parameter(zsolar1=0.517) +! parameter(zsolar2=0.483) ! !* 2.2) namelist parameters (see *sub* radini) ! + real :: starbbtemp = 5607.9 ! Star's blackbody surface temperature (K) + real :: gsol0 = 1367.0 ! solar constant (set in planet module) real :: solclat = 1.0 ! cos of lat of insolation if ncstsol=1 real :: solcdec = 1.0 ! cos of dec of insolation if ncstsol=1 @@ -55,6 +59,8 @@ integer :: nfixed = 0 ! Switch for fixed zenith angle (0/1=no/yes) real :: fixedlon = 0.0 ! Longitude of fixed solar zenith real :: slowdown = 1.0 ! Factor by which to change diurnal insolation cycle + + integer :: nstartemp = 0 ! Switch for using the star's bb temp to determine sw (0/1) real :: rcl1(3)=(/0.15,0.30,0.60/) ! cloud albedos spectral range 1 real :: rcl2(3)=(/0.15,0.30,0.60/) ! cloud albedos spectral range 2 @@ -68,9 +74,13 @@ real :: gmu1(NHOR) ! cosine of solar zenit angle real :: dqo3(NHOR,NLEV) = 0.0 ! ozon concentration (kg/kg) real :: dqco2(NHOR,NLEV) = 0.0 ! co2 concentration (ppmv) + real :: dtdtlwr(NHOR,NLEV) ! lwr temperature tendencies + real :: dtdtswr(NHOR,NLEV) ! swr temperature tendencies real, allocatable :: dqo3cl(:,:,:) ! climatological O3 (used if NO3=2) + real :: zsolars(2) = 0.0 ! Container for storing solar constants + ! !* 2.4) scalars ! @@ -79,6 +89,9 @@ real :: time4rad = 0. ! CPU time for radiation real :: time4swr = 0. ! CPU time for short wave radiation real :: time4lwr = 0. ! CPU time for long wave radiation + + real :: zsolar1 = 0.517 ! spectral partitioning 1 (wl < 0.75mue) + real :: zsolar2 = 0.483 ! spectral partitioning 2 (wl > 0.75mue) ! ! 2.5 orbital parameters @@ -116,6 +129,567 @@ ! radiation subroutines ! +! =================== +! SUBROUTINE SOLARINI +! =================== + +! +! Spectral snow reflectance sourced from the ECOSTRESS spectral library (JPL,JHU,USGS): +! +! Meerdink, S. K., Hook, S. J., Abbott, E.A., & Roberts, D.A. (in prep). +! The ECOSTRESS Spectral Library 1.0. +! Baldridge, A. M., S.J. Hook, C.I. Grove and G. Rivera, 2009.. +! The ASTER Spectral Library Version 2.0. Remote Sensing of Environment, vol 113, pp. 711-715 +! + subroutine solarini + use radmod + +! parameter(planckh = 6.62607004e-34) +! parameter(boltzk = 1.38064852e-23 ) +! parameter(cc = 299792458.0 ) + parameter(const = 0.0143877735383) !hc/k + + real :: wv1(1024) !Wavelengths in meters up to 0.75 microns + real :: wv2(1024) !Wavelength in meters starting at 0.75 microns + real :: wvm1(1024) !Wavelengths in microns up to 0.75 microns + real :: wvm2(1024) !Wavelength in microns starting at 0.75 microns + real :: bb1(1024) !Planck function for x<0.75 microns + real :: bb2(1024) !Planck function for x>0.75 microns + real :: bb3(965) !Planck function for albedo wavelengths + + real dl1,dl2,hinge,const1,const2,z1,z2,znet,wmin,lwmin + integer k + + real :: snowalbedos(965) = 0.0 !Combined snow spectral albedos + + real :: wavelengths(965) = (/0.3400,0.3500,0.3600,0.3700,0.3800,0.3900,0.4000,0.4100,0.4200, & + & 0.4300,0.4400,0.4500,0.4600,0.4700,0.4800,0.4900,0.5000,0.5100, & + & 0.5200,0.5300,0.5400,0.5500,0.5600,0.5700,0.5800,0.5900,0.6000, & + & 0.6100,0.6200,0.6300,0.6400,0.6500,0.6600,0.6700,0.6800,0.6900, & + & 0.7000,0.7100,0.7200,0.7300,0.7400,0.7500,0.7600,0.7700,0.7800, & + & 0.7900,0.8000,0.8200,0.8400,0.8600,0.8800,0.9000,0.9200,0.9400, & + & 0.9600,0.9800,1.0000,1.0200,1.0400,1.0600,1.0800,1.1000,1.1200, & + & 1.1400,1.1600,1.1800,1.2000,1.2200,1.2400,1.2600,1.2800,1.3000, & + & 1.3200,1.3400,1.3600,1.3800,1.4000,1.4200,1.4400,1.4600,1.4800, & + & 1.5000,1.5200,1.5400,1.5600,1.5800,1.6000,1.6200,1.6400,1.6600, & + & 1.6800,1.7000,1.7200,1.7400,1.7600,1.7800,1.8000,1.8200,1.8400, & + & 1.8600,1.8800,1.9000,1.9200,1.9400,1.9600,1.9800,2.0000,2.0200, & + & 2.0400,2.0600,2.0800,2.1000,2.1200,2.1400,2.1600,2.1800,2.2000, & + & 2.2200,2.2400,2.2500,2.2600,2.2700,2.2800,2.2900,2.3000,2.3100, & + & 2.3200,2.3300,2.3400,2.3500,2.3600,2.3700,2.3800,2.3900,2.4000, & + & 2.4100,2.4200,2.4300,2.4400,2.4500,2.4600,2.4700,2.4800,2.4900, & + & 2.5000,2.5100,2.5200,2.5300,2.5400,2.5500,2.5600,2.5700,2.5800, & + & 2.5900,2.6000,2.6100,2.6200,2.6300,2.6400,2.6500,2.6600,2.6700, & + & 2.6800,2.6900,2.7000,2.7100,2.7200,2.7300,2.7400,2.7500,2.7600, & + & 2.7700,2.7800,2.7900,2.8000,2.8100,2.8200,2.8300,2.8400,2.8500, & + & 2.8600,2.8700,2.8800,2.8900,2.9000,2.9100,2.9200,2.9300,2.9400, & + & 2.9500,2.9600,2.9700,2.9800,2.9900,3.0000,3.0100,3.0200,3.0300, & + & 3.0400,3.0500,3.0600,3.0700,3.0800,3.0900,3.1000,3.1100,3.1200, & + & 3.1300,3.1400,3.1500,3.1600,3.1700,3.1800,3.1900,3.2000,3.2100, & + & 3.2200,3.2300,3.2400,3.2500,3.2600,3.2700,3.2800,3.2900,3.3000, & + & 3.3100,3.3200,3.3300,3.3400,3.3500,3.3600,3.3700,3.3800,3.3900, & + & 3.4000,3.4100,3.4200,3.4300,3.4400,3.4500,3.4600,3.4700,3.4800, & + & 3.4900,3.5000,3.5100,3.5200,3.5300,3.5400,3.5500,3.5600,3.5700, & + & 3.5800,3.5900,3.6000,3.6100,3.6200,3.6300,3.6400,3.6500,3.6600, & + & 3.6700,3.6800,3.6900,3.7000,3.7100,3.7200,3.7300,3.7400,3.7500, & + & 3.7600,3.7700,3.7800,3.7900,3.8000,3.8100,3.8200,3.8300,3.8400, & + & 3.8500,3.8600,3.8700,3.8800,3.8900,3.9000,3.9100,3.9200,3.9300, & + & 3.9400,3.9500,3.9600,3.9700,3.9800,3.9900,4.0000,4.0100,4.0200, & + & 4.0300,4.0400,4.0500,4.0600,4.0700,4.0800,4.0900,4.1000,4.1100, & + & 4.1200,4.1300,4.1400,4.1500,4.1600,4.1700,4.1800,4.1900,4.2000, & + & 4.2100,4.2200,4.2300,4.2400,4.2500,4.2600,4.2700,4.2800,4.2900, & + & 4.3000,4.3100,4.3200,4.3300,4.3400,4.3500,4.3600,4.3700,4.3800, & + & 4.3900,4.4000,4.4100,4.4200,4.4300,4.4400,4.4500,4.4600,4.4700, & + & 4.4800,4.4900,4.5000,4.5100,4.5200,4.5300,4.5400,4.5500,4.5600, & + & 4.5700,4.5800,4.5900,4.6000,4.6100,4.6200,4.6300,4.6400,4.6500, & + & 4.6600,4.6700,4.6800,4.6900,4.7000,4.7100,4.7200,4.7300,4.7400, & + & 4.7500,4.7600,4.7700,4.7800,4.7900,4.8000,4.8100,4.8200,4.8300, & + & 4.8400,4.8500,4.8600,4.8700,4.8800,4.8900,4.9000,4.9100,4.9200, & + & 4.9300,4.9400,4.9500,4.9600,4.9700,4.9800,4.9900,5.0000,5.0100, & + & 5.0200,5.0300,5.0400,5.0500,5.0600,5.0700,5.0800,5.0900,5.1000, & + & 5.1100,5.1200,5.1300,5.1400,5.1500,5.1600,5.1700,5.1800,5.1900, & + & 5.2000,5.2100,5.2200,5.2300,5.2400,5.2500,5.2600,5.2700,5.2800, & + & 5.2900,5.3000,5.3100,5.3200,5.3300,5.3400,5.3500,5.3600,5.3700, & + & 5.3800,5.3900,5.4000,5.4100,5.4200,5.4300,5.4400,5.4500,5.4600, & + & 5.4700,5.4800,5.4900,5.5000,5.5100,5.5200,5.5300,5.5400,5.5500, & + & 5.5600,5.5700,5.5800,5.5900,5.6000,5.6100,5.6200,5.6300,5.6400, & + & 5.6500,5.6600,5.6700,5.6800,5.6900,5.7000,5.7100,5.7200,5.7300, & + & 5.7400,5.7500,5.7600,5.7700,5.7800,5.7900,5.8000,5.8100,5.8200, & + & 5.8300,5.8400,5.8500,5.8600,5.8700,5.8800,5.8900,5.9000,5.9100, & + & 5.9200,5.9300,5.9400,5.9500,5.9600,5.9700,5.9800,5.9900,6.0000, & + & 6.0100,6.0200,6.0300,6.0400,6.0500,6.0600,6.0700,6.0800,6.0900, & + & 6.1000,6.1100,6.1200,6.1300,6.1400,6.1500,6.1600,6.1700,6.1800, & + & 6.1900,6.2000,6.2100,6.2200,6.2300,6.2400,6.2500,6.2600,6.2700, & + & 6.2800,6.2900,6.3000,6.3100,6.3200,6.3300,6.3400,6.3500,6.3600, & + & 6.3700,6.3800,6.3900,6.4000,6.4100,6.4200,6.4300,6.4400,6.4500, & + & 6.4600,6.4700,6.4800,6.4900,6.5000,6.5100,6.5200,6.5300,6.5400, & + & 6.5500,6.5600,6.5700,6.5800,6.5900,6.6000,6.6100,6.6200,6.6300, & + & 6.6400,6.6500,6.6600,6.6700,6.6800,6.6900,6.7000,6.7100,6.7200, & + & 6.7300,6.7400,6.7500,6.7600,6.7700,6.7800,6.7900,6.8000,6.8100, & + & 6.8200,6.8300,6.8400,6.8500,6.8600,6.8700,6.8800,6.8900,6.9000, & + & 6.9100,6.9200,6.9300,6.9400,6.9500,6.9600,6.9700,6.9800,6.9900, & + & 7.0000,7.0100,7.0200,7.0300,7.0400,7.0500,7.0600,7.0700,7.0800, & + & 7.0900,7.1000,7.1100,7.1200,7.1300,7.1400,7.1500,7.1600,7.1700, & + & 7.1800,7.1900,7.2000,7.2100,7.2200,7.2300,7.2400,7.2500,7.2600, & + & 7.2700,7.2800,7.2900,7.3000,7.3100,7.3200,7.3300,7.3400,7.3500, & + & 7.3600,7.3700,7.3800,7.4000,7.4100,7.4200,7.4300,7.4400,7.4500, & + & 7.4600,7.4700,7.4800,7.4900,7.5000,7.5100,7.5200,7.5400,7.5500, & + & 7.5600,7.5700,7.5800,7.5900,7.6000,7.6100,7.6200,7.6300,7.6500, & + & 7.6600,7.6700,7.6800,7.6900,7.7000,7.7100,7.7300,7.7400,7.7500, & + & 7.7600,7.7700,7.7800,7.8000,7.8100,7.8200,7.8300,7.8400,7.8500, & + & 7.8700,7.8800,7.8900,7.9000,7.9100,7.9300,7.9400,7.9500,7.9600, & + & 7.9800,7.9900,8.0000,8.0100,8.0300,8.0400,8.0500,8.0600,8.0800, & + & 8.0900,8.1000,8.1100,8.1300,8.1400,8.1500,8.1600,8.1800,8.1900, & + & 8.2000,8.2200,8.2300,8.2400,8.2600,8.2700,8.2800,8.2900,8.3100, & + & 8.3200,8.3300,8.3500,8.3600,8.3800,8.3900,8.4000,8.4200,8.4300, & + & 8.4400,8.4600,8.4700,8.4800,8.5000,8.5100,8.5300,8.5400,8.5500, & + & 8.5700,8.5800,8.6000,8.6100,8.6300,8.6400,8.6500,8.6700,8.6800, & + & 8.7000,8.7100,8.7300,8.7400,8.7600,8.7700,8.7900,8.8000,8.8200, & + & 8.8300,8.8500,8.8600,8.8800,8.8900,8.9100,8.9200,8.9400,8.9500, & + & 8.9700,8.9800,9.0000,9.0200,9.0300,9.0500,9.0600,9.0800,9.1000, & + & 9.1100,9.1300,9.1400,9.1600,9.1800,9.1900,9.2100,9.2200,9.2400, & + & 9.2600,9.2700,9.2900,9.3100,9.3200,9.3400,9.3600,9.3700,9.3900, & + & 9.4100,9.4300,9.4400,9.4600,9.4800,9.4900,9.5100,9.5300,9.5500, & + & 9.5600,9.5800,9.6000,9.6200,9.6400,9.6500,9.6700,9.6900,9.7100, & + & 9.7300,9.7400,9.7600,9.7800,9.8000,9.8200,9.8400,9.8600,9.8700, & + & 9.8900,9.9100,9.9300,9.9500,9.9700,9.9900,10.0100,10.0300,10.0500, & + & 10.0700,10.0900,10.1100,10.1300,10.1500,10.1700,10.1800,10.2100,10.2300, & + & 10.2500,10.2700,10.2900,10.3100,10.3300,10.3500,10.3700,10.3900,10.4100, & + & 10.4300,10.4500,10.4700,10.4900,10.5200,10.5400,10.5600,10.5800,10.6000, & + & 10.6200,10.6500,10.6700,10.6900,10.7100,10.7300,10.7600,10.7800,10.8000, & + & 10.8200,10.8500,10.8700,10.8900,10.9100,10.9400,10.9600,10.9800,11.0100, & + & 11.0300,11.0500,11.0800,11.1000,11.1200,11.1500,11.1700,11.2000,11.2200, & + & 11.2500,11.2700,11.2900,11.3200,11.3400,11.3700,11.3900,11.4200,11.4400, & + & 11.4700,11.4900,11.5200,11.5500,11.5700,11.6000,11.6200,11.6500,11.6800, & + & 11.7000,11.7300,11.7600,11.7800,11.8100,11.8400,11.8600,11.8900,11.9200, & + & 11.9500,11.9700,12.0000,12.0300,12.0600,12.0800,12.1100,12.1400,12.1700, & + & 12.2000,12.2300,12.2600,12.2800,12.3100,12.3400,12.3700,12.4000,12.4300, & + & 12.4600,12.4900,12.5200,12.5500,12.5800,12.6100,12.6400,12.6800,12.7100, & + & 12.7400,12.7700,12.8000,12.8300,12.8600,12.9000,12.9300,12.9600,12.9900, & + & 13.0300,13.0600,13.0900,13.1200,13.1600,13.1900,13.2200,13.2600,13.2900, & + & 13.3300,13.3600,13.4000,13.4300,13.4700,13.5000,13.5400,13.5700,13.6100, & + & 13.6400,13.6800,13.7100,13.7500,13.7900,13.8200,13.8600,13.9000,13.9400, & + & 13.9700,14.0100/) + + real :: csnowalb(965) = (/96.0980,96.2900,96.4230,96.5490,96.7270,96.8430,96.9870,97.1520,97.2650, & + & 97.4450,97.5050,97.8490,97.9080,97.8730,97.7810,97.7360,97.6340,97.5820, & + & 97.5280,97.4620,97.2960,97.1630,97.1230,97.0590,96.8980,96.7230,96.3570, & + & 95.5080,95.6910,95.2990,94.9000,94.5670,94.2170,93.6110,93.6190,93.1990, & + & 92.5050,92.0360,91.5250,91.2130,90.6720,89.9820,89.1350,88.1980,87.0520, & + & 86.4820,85.8460,85.2430,85.2360,82.8540,79.3290,77.3930,76.2460,74.9770, & + & 71.6990,67.0850,60.8970,57.4940,57.1240,60.0210,62.2260,62.7290,61.9590, & + & 59.0470,54.5580,47.3760,41.5390,34.0440,31.5480,30.3300,30.5110,30.8370, & + & 31.1360,31.2110,30.2540,28.7190,25.1350,11.2590,2.3120,1.0520,0.6190, & + & 0.5090,0.5610,0.6340,0.8010,1.0100,1.1870,1.4250,1.8840,2.4420, & + & 3.0950,3.7230,4.2260,4.7940,5.3750,6.1670,7.0920,8.1850,9.2150, & + & 7.5990,3.0310,1.2170,0.5350,0.3380,0.2420,0.2040,0.2000,0.2070, & + & 0.2260,0.2610,0.3140,0.4210,0.6530,1.0260,1.4710,1.9870,2.6210, & + & 3.3170,3.6960,3.6770,3.5800,3.3610,3.0840,2.7980,2.5400,2.3020, & + & 2.1370,2.0200,1.9140,1.8110,1.7500,1.7000,1.6820,1.6370,1.6010, & + & 1.6030,1.5740,1.5560,1.5080,1.4980,1.4850,1.4660,1.4530,1.4370, & + & 1.3800,1.3590,1.3670,1.3510,1.3070,1.3150,1.3290,1.2940,1.3010, & + & 1.2920,1.2700,1.2390,1.2150,1.1570,1.1370,1.0880,1.0650,1.0010, & + & 0.9250,0.8920,0.8090,0.7620,0.7060,0.6460,0.6050,0.5840,0.5840, & + & 0.5250,0.4850,0.4820,0.4300,0.4150,0.4000,0.3990,0.3940,0.4130, & + & 0.4420,0.4660,0.5250,0.6480,0.7990,0.9970,1.1950,1.4360,1.7020, & + & 1.9030,2.0990,2.2790,2.4280,2.6390,2.8780,3.2290,3.6490,4.1500, & + & 4.6890,5.3580,5.9130,6.4650,6.7620,6.8220,6.6880,6.3590,6.0790, & + & 5.8940,5.7770,5.7110,5.6730,5.6280,5.5760,5.4410,5.3060,5.0970, & + & 4.8930,4.6840,4.4740,4.3040,4.1440,3.9930,3.8770,3.7490,3.6490, & + & 3.5630,3.4740,3.3850,3.3020,3.2380,3.2010,3.1450,3.0470,2.9840, & + & 2.9590,2.9320,2.8880,2.8330,2.7710,2.7260,2.7030,2.6570,2.6100, & + & 2.5860,2.5670,2.5300,2.4780,2.4840,2.4590,2.4290,2.4070,2.3850, & + & 2.3810,2.3670,2.3480,2.3140,2.2960,2.3010,2.2700,2.2490,2.2400, & + & 2.2150,2.2070,2.1990,2.1840,2.1490,2.1610,2.1580,2.1240,2.1190, & + & 2.1150,2.0950,2.0830,2.0800,2.0720,2.0540,2.0460,2.0490,2.0530, & + & 2.0370,2.0050,2.0050,2.0090,1.9850,1.9670,1.9620,1.9640,1.9610, & + & 1.9600,1.9560,1.9350,1.9280,1.9290,1.9230,1.9240,1.9160,1.9020, & + & 1.9000,1.8830,1.8570,1.8670,1.8880,1.8740,1.8520,1.8570,1.8600, & + & 1.8480,1.8440,1.8440,1.8340,1.8330,1.8320,1.8230,1.8330,1.8350, & + & 1.8120,1.7980,1.8020,1.8090,1.8040,1.8010,1.7990,1.7940,1.7950, & + & 1.8000,1.7930,1.7810,1.7830,1.7940,1.7940,1.7850,1.7820,1.7790, & + & 1.7830,1.7910,1.7960,1.7980,1.7830,1.7810,1.8070,1.8270,1.8230, & + & 1.8150,1.8170,1.8190,1.8160,1.8190,1.8280,1.8320,1.8340,1.8400, & + & 1.8420,1.8410,1.8420,1.8430,1.8510,1.8640,1.8700,1.8700,1.8640, & + & 1.8570,1.8520,1.8490,1.8490,1.8510,1.8480,1.8370,1.8250,1.8230, & + & 1.8230,1.8150,1.8110,1.8190,1.8290,1.8220,1.8090,1.8030,1.8010, & + & 1.7990,1.7970,1.7950,1.7920,1.7880,1.7870,1.7880,1.7850,1.7800, & + & 1.7810,1.7860,1.7840,1.7740,1.7640,1.7540,1.7450,1.7420,1.7430, & + & 1.7470,1.7490,1.7440,1.7360,1.7330,1.7330,1.7300,1.7210,1.7080, & + & 1.7010,1.7060,1.7150,1.7100,1.6950,1.6830,1.6830,1.6840,1.6770, & + & 1.6630,1.6550,1.6510,1.6510,1.6530,1.6560,1.6580,1.6530,1.6460, & + & 1.6440,1.6440,1.6400,1.6270,1.6190,1.6260,1.6360,1.6430,1.6390, & + & 1.6260,1.6110,1.6010,1.5920,1.5860,1.5840,1.5800,1.5730,1.5650, & + & 1.5600,1.5570,1.5580,1.5670,1.5740,1.5710,1.5530,1.5360,1.5320, & + & 1.5360,1.5390,1.5360,1.5280,1.5170,1.5120,1.5140,1.5180,1.5150, & + & 1.5040,1.4940,1.4920,1.4940,1.4890,1.4740,1.4590,1.4550,1.4630, & + & 1.4700,1.4680,1.4620,1.4600,1.4640,1.4730,1.4820,1.4830,1.4770, & + & 1.4720,1.4710,1.4730,1.4700,1.4620,1.4490,1.4370,1.4350,1.4400, & + & 1.4520,1.4610,1.4640,1.4660,1.4710,1.4820,1.4940,1.4990,1.4970, & + & 1.4940,1.4940,1.4950,1.4970,1.4980,1.5050,1.5200,1.5370,1.5500, & + & 1.5560,1.5570,1.5550,1.5540,1.5510,1.5490,1.5500,1.5540,1.5590, & + & 1.5660,1.5710,1.5740,1.5750,1.5740,1.5710,1.5700,1.5750,1.5830, & + & 1.5970,1.6070,1.6080,1.6050,1.6020,1.6020,1.6050,1.6130,1.6200, & + & 1.6220,1.6170,1.6060,1.5980,1.5940,1.6020,1.6240,1.6410,1.6550, & + & 1.6640,1.6640,1.6630,1.6640,1.6650,1.6650,1.6610,1.6570,1.6540, & + & 1.6540,1.6610,1.6690,1.6770,1.6820,1.6840,1.6810,1.6780,1.6770, & + & 1.6780,1.6820,1.6870,1.6900,1.6920,1.6950,1.6990,1.7030,1.7040, & + & 1.6990,1.6930,1.6870,1.6820,1.6780,1.6740,1.6690,1.6640,1.6600, & + & 1.6580,1.6570,1.6600,1.6650,1.6730,1.6810,1.6900,1.6920,1.6900, & + & 1.6870,1.6830,1.6810,1.6800,1.6800,1.6820,1.6860,1.6890,1.6930, & + & 1.6950,1.6950,1.6930,1.6900,1.6890,1.6890,1.6890,1.6870,1.6820, & + & 1.6760,1.6690,1.6640,1.6620,1.6640,1.6680,1.6730,1.6750,1.6740, & + & 1.6710,1.6660,1.6610,1.6580,1.6550,1.6520,1.6490,1.6470,1.6470, & + & 1.6470,1.6460,1.6430,1.6390,1.6330,1.6260,1.6210,1.6170,1.6150, & + & 1.6140,1.6140,1.6160,1.6190,1.6250,1.6290,1.6300,1.6250,1.6140, & + & 1.5980,1.5800,1.5670,1.5620,1.5660,1.5770,1.5900,1.5990,1.6010, & + & 1.5960,1.5850,1.5730,1.5620,1.5550,1.5520,1.5530,1.5550,1.5580, & + & 1.5590,1.5570,1.5530,1.5490,1.5460,1.5440,1.5450,1.5460,1.5470, & + & 1.5470,1.5460,1.5460,1.5450,1.5450,1.5440,1.5430,1.5380,1.5330, & + & 1.5260,1.5220,1.5200,1.5230,1.5280,1.5320,1.5340,1.5330,1.5310, & + & 1.5290,1.5300,1.5330,1.5360,1.5390,1.5390,1.5380,1.5360,1.5330, & + & 1.5290,1.5260,1.5240,1.5250,1.5300,1.5370,1.5450,1.5510,1.5530, & + & 1.5500,1.5450,1.5390,1.5320,1.5250,1.5180,1.5120,1.5090,1.5100, & + & 1.5140,1.5190,1.5240,1.5260,1.5260,1.5250,1.5230,1.5240,1.5260, & + & 1.5280,1.5270,1.5210,1.5110,1.4970,1.4830,1.4720,1.4660,1.4650, & + & 1.4680,1.4740,1.4800,1.4850,1.4880,1.4900,1.4900,1.4900,1.4900, & + & 1.4880,1.4850,1.4790,1.4720,1.4660,1.4620,1.4620,1.4650,1.4710, & + & 1.4750,1.4730,1.4660,1.4560,1.4480,1.4460,1.4510,1.4620,1.4700, & + & 1.4700,1.4610,1.4450,1.4260,1.4110,1.4050,1.4060,1.4120,1.4180, & + & 1.4230,1.4230,1.4210,1.4180,1.4140,1.4090,1.4040,1.3970,1.3900, & + & 1.3840,1.3790,1.3740,1.3690,1.3610,1.3510,1.3420,1.3330,1.3270, & + & 1.3230,1.3210,1.3160,1.3080,1.2980,1.2830,1.2660,1.2470,1.2300, & + & 1.2160,1.2090,1.2070,1.2110,1.2160,1.2170,1.2110,1.2000,1.1820, & + & 1.1610,1.1390,1.1170,1.0970,1.0790,1.0630,1.0490,1.0380,1.0280, & + & 1.0190,1.0080,0.9950,0.9810,0.9660,0.9530,0.9430,0.9340,0.9280, & + & 0.9220,0.9190,0.9150,0.9090,0.9000,0.8860,0.8660,0.8400,0.8130, & + & 0.7900,0.7760,0.7710,0.7730,0.7780,0.7790,0.7760,0.7680,0.7580, & + & 0.7520,0.7520,0.7560,0.7630,0.7730,0.7870,0.8020,0.8210,0.8450, & + & 0.8730,0.9050,0.9410,0.9840,1.0310,1.0810,1.1290,1.1760,1.2210, & + & 1.2650,1.3120,1.3630,1.4180,1.4750,1.5310,1.5880,1.6460,1.7070, & + & 1.7710,1.8370,1.9030,1.9660,2.0250,2.0820,2.1390,2.1980,2.2590, & + & 2.3200,2.3810,2.4400,2.4950,2.5490,2.6030,2.6580,2.7140,2.7720, & + & 2.8310,2.8910,2.9500,3.0030,3.0490,3.0880,3.1210,3.1510,3.1840, & + & 3.2200,3.2600,3.3040,3.3480,3.3930,3.4370,3.4820,3.5280,3.5770, & + & 3.6290,3.6840,3.7400,3.7950,3.8480,3.8950,3.9400,3.9830,4.0240, & + & 4.0640,4.1060,4.1460,4.1840,4.2180,4.2480,4.2750,4.2960,4.3130, & + & 4.3300,4.3470,4.3670,4.3890,4.4130,4.4390,4.4630,4.4850,4.5050, & + & 4.5230,4.5360,4.5420,4.5400,4.5330,4.5250,4.5190,4.5170,4.5200, & + & 4.5230,4.5240,4.5230,4.5190,4.5110,4.5000,4.4850,4.4660,4.4440, & + & 4.4260,4.4210,4.4280,4.4450,4.4640,4.4760,4.4770,4.4660,4.4490, & + & 4.4310,4.4120,4.3870,4.3540,4.3120,4.2670,4.2260,4.1950,4.1810, & + & 4.1740,4.1600 /) + + real :: msnowalb(965) = (/97.3810,97.4320,97.5420,97.6360,97.8130,97.8360,97.9740,98.0590,98.2080, & + & 98.2900,98.3790,98.5380,98.5360,98.5780,98.5650,98.5220,98.4490,98.3590, & + & 98.3240,98.2340,98.1110,98.1510,98.0970,98.0130,97.9010,97.6920,97.5250, & + & 97.2970,97.0760,96.8090,96.5850,95.9750,96.0250,95.7930,95.6150,95.3290, & + & 94.9490,94.5520,94.1700,93.9720,93.6170,93.0600,92.4990,91.8540,91.1910, & + & 90.3890,90.0660,89.9020,89.7670,88.0740,85.4770,83.9870,83.3810,82.3000, & + & 79.7740,76.1730,72.4830,68.6590,67.9870,70.8230,72.4330,72.8880,72.3670, & + & 70.0170,66.2320,60.1390,55.0300,48.2160,44.6580,44.5010,44.3870,44.9400, & + & 45.1580,45.0200,43.8600,42.6350,38.9960,21.5660,6.4470,3.2650,1.9890, & + & 1.6340,1.8050,2.0320,2.5540,3.1620,3.6430,4.2570,5.4520,6.7560, & + & 8.1060,9.4780,10.5090,11.7060,12.7990,13.9730,15.5850,17.2010,19.1400, & + & 16.5080,8.1560,3.7400,1.7400,1.0750,0.7150,0.5580,0.5410,0.5740, & + & 0.6560,0.7920,0.9950,1.3730,2.1160,3.1900,4.4140,5.6970,7.0460, & + & 8.6660,9.6390,9.5940,9.1440,8.3630,7.3170,6.2050,5.2290,4.4140, & + & 3.7100,3.2060,2.8250,2.5330,2.3520,2.2160,2.1170,2.0610,1.9720, & + & 1.9510,1.8750,1.8350,1.8090,1.7320,1.6840,1.6250,1.5760,1.5420, & + & 1.5030,1.4790,1.4690,1.4600,1.4210,1.4290,1.4430,1.4510,1.4390, & + & 1.4190,1.3970,1.3610,1.3150,1.2530,1.1880,1.1020,1.0230,0.9270, & + & 0.8490,0.7560,0.6740,0.5800,0.5480,0.5070,0.4470,0.4310,0.4230, & + & 0.3880,0.3780,0.3570,0.3310,0.3130,0.3040,0.2910,0.2740,0.2760, & + & 0.2910,0.3120,0.3660,0.4360,0.5410,0.6910,0.8110,0.9300,1.0700, & + & 1.1630,1.2620,1.3780,1.4420,1.5570,1.6740,1.8370,2.0640,2.3210, & + & 2.6550,3.0310,3.3710,3.6570,3.8060,3.7920,3.6440,3.4540,3.2730, & + & 3.1290,3.0500,3.0170,3.0000,2.9650,2.9180,2.8580,2.7630,2.6530, & + & 2.5330,2.4230,2.3380,2.2430,2.1580,2.0780,2.0320,1.9920,1.9490, & + & 1.9060,1.8750,1.8390,1.8180,1.8090,1.7800,1.7620,1.7490,1.7380, & + & 1.7240,1.6960,1.7000,1.7070,1.6810,1.6570,1.6590,1.6550,1.6480, & + & 1.6400,1.6310,1.6300,1.6200,1.6030,1.5890,1.5840,1.5840,1.5760, & + & 1.5690,1.5620,1.5500,1.5430,1.5420,1.5310,1.5110,1.5070,1.5110, & + & 1.5060,1.4940,1.4830,1.4660,1.4530,1.4590,1.4600,1.4450,1.4280, & + & 1.4170,1.4070,1.4000,1.3970,1.3890,1.3800,1.3740,1.3600,1.3470, & + & 1.3460,1.3450,1.3320,1.3170,1.3120,1.3110,1.3060,1.2990,1.2950, & + & 1.2900,1.2770,1.2670,1.2640,1.2550,1.2460,1.2460,1.2460,1.2390, & + & 1.2280,1.2160,1.2100,1.2100,1.2120,1.2090,1.2020,1.1930,1.1900, & + & 1.1910,1.1920,1.1920,1.1880,1.1860,1.1840,1.1750,1.1610,1.1510, & + & 1.1460,1.1450,1.1450,1.1460,1.1460,1.1440,1.1420,1.1410,1.1400, & + & 1.1380,1.1330,1.1280,1.1210,1.1150,1.1120,1.1120,1.1100,1.1080, & + & 1.1070,1.1070,1.1060,1.1080,1.1190,1.1280,1.1270,1.1200,1.1190, & + & 1.1250,1.1290,1.1280,1.1280,1.1320,1.1370,1.1430,1.1500,1.1570, & + & 1.1620,1.1660,1.1680,1.1700,1.1740,1.1800,1.1830,1.1830,1.1840, & + & 1.1840,1.1840,1.1810,1.1780,1.1760,1.1790,1.1850,1.1920,1.1950, & + & 1.1950,1.1910,1.1860,1.1830,1.1820,1.1830,1.1810,1.1770,1.1730, & + & 1.1710,1.1750,1.1830,1.1900,1.1930,1.1900,1.1820,1.1780,1.1780, & + & 1.1810,1.1830,1.1850,1.1870,1.1860,1.1810,1.1750,1.1720,1.1700, & + & 1.1710,1.1710,1.1710,1.1700,1.1670,1.1650,1.1620,1.1570,1.1520, & + & 1.1470,1.1450,1.1420,1.1370,1.1290,1.1210,1.1180,1.1170,1.1170, & + & 1.1140,1.1110,1.1070,1.1040,1.1020,1.1020,1.1020,1.0990,1.0950, & + & 1.0900,1.0850,1.0800,1.0760,1.0740,1.0730,1.0720,1.0720,1.0730, & + & 1.0720,1.0690,1.0640,1.0580,1.0520,1.0490,1.0490,1.0500,1.0480, & + & 1.0460,1.0430,1.0390,1.0320,1.0260,1.0200,1.0130,1.0080,1.0050, & + & 1.0040,1.0040,1.0050,1.0080,1.0100,1.0080,1.0020,0.9900,0.9820, & + & 0.9750,0.9700,0.9670,0.9630,0.9600,0.9560,0.9510,0.9450,0.9380, & + & 0.9340,0.9330,0.9340,0.9360,0.9360,0.9340,0.9300,0.9260,0.9220, & + & 0.9220,0.9260,0.9310,0.9350,0.9370,0.9340,0.9280,0.9190,0.9100, & + & 0.9040,0.8990,0.8970,0.8960,0.8960,0.8960,0.8980,0.9020,0.9060, & + & 0.9090,0.9090,0.9070,0.9000,0.8910,0.8850,0.8820,0.8840,0.8890, & + & 0.8930,0.8970,0.9000,0.9010,0.9030,0.9040,0.9050,0.9070,0.9070, & + & 0.9080,0.9090,0.9090,0.9110,0.9120,0.9130,0.9140,0.9150,0.9160, & + & 0.9170,0.9180,0.9190,0.9200,0.9220,0.9230,0.9240,0.9270,0.9280, & + & 0.9290,0.9300,0.9320,0.9340,0.9380,0.9460,0.9540,0.9590,0.9640, & + & 0.9690,0.9730,0.9750,0.9770,0.9800,0.9830,0.9850,0.9870,0.9890, & + & 0.9930,0.9950,0.9960,0.9970,0.9980,0.9970,0.9970,0.9960,0.9960, & + & 0.9960,0.9960,0.9970,0.9980,0.9990,1.0000,1.0010,1.0020,1.0010, & + & 1.0000,0.9980,0.9970,0.9950,0.9940,0.9930,0.9920,0.9900,0.9890, & + & 0.9880,0.9860,0.9850,0.9840,0.9850,0.9870,0.9930,0.9990,1.0030, & + & 1.0060,1.0090,1.0100,1.0100,1.0100,1.0100,1.0110,1.0130,1.0140, & + & 1.0170,1.0190,1.0200,1.0200,1.0180,1.0160,1.0130,1.0110,1.0110, & + & 1.0120,1.0150,1.0190,1.0230,1.0270,1.0300,1.0300,1.0300,1.0270, & + & 1.0240,1.0210,1.0170,1.0140,1.0120,1.0100,1.0080,1.0060,1.0030, & + & 1.0010,0.9980,0.9940,0.9910,0.9880,0.9860,0.9840,0.9820,0.9800, & + & 0.9790,0.9770,0.9760,0.9750,0.9730,0.9720,0.9710,0.9700,0.9690, & + & 0.9690,0.9680,0.9680,0.9670,0.9670,0.9660,0.9650,0.9630,0.9600, & + & 0.9580,0.9550,0.9530,0.9520,0.9510,0.9500,0.9500,0.9500,0.9500, & + & 0.9500,0.9510,0.9510,0.9520,0.9520,0.9520,0.9520,0.9510,0.9490, & + & 0.9470,0.9460,0.9440,0.9440,0.9440,0.9450,0.9470,0.9490,0.9510, & + & 0.9530,0.9540,0.9540,0.9540,0.9530,0.9510,0.9490,0.9480,0.9460, & + & 0.9440,0.9430,0.9410,0.9410,0.9410,0.9410,0.9410,0.9410,0.9400, & + & 0.9390,0.9380,0.9370,0.9360,0.9360,0.9370,0.9390,0.9420,0.9460, & + & 0.9490,0.9530,0.9580,0.9610,0.9630,0.9640,0.9630,0.9600,0.9560, & + & 0.9500,0.9430,0.9350,0.9280,0.9220,0.9180,0.9160,0.9160,0.9170, & + & 0.9200,0.9230,0.9260,0.9270,0.9260,0.9230,0.9190,0.9140,0.9100, & + & 0.9060,0.9020,0.9010,0.8990,0.8990,0.8980,0.8980,0.8990,0.9000, & + & 0.9010,0.9010,0.9020,0.9030,0.9020,0.9010,0.8990,0.8970,0.8950, & + & 0.8940,0.8930,0.8920,0.8910,0.8900,0.8890,0.8880,0.8870,0.8870, & + & 0.8860,0.8860,0.8840,0.8820,0.8770,0.8710,0.8630,0.8540,0.8460, & + & 0.8400,0.8350,0.8330,0.8350,0.8380,0.8430,0.8500,0.8570,0.8630, & + & 0.8670,0.8690,0.8680,0.8640,0.8580,0.8500,0.8410,0.8310,0.8210, & + & 0.8110,0.8030,0.7960,0.7900,0.7840,0.7790,0.7750,0.7700,0.7650, & + & 0.7600,0.7530,0.7460,0.7370,0.7260,0.7160,0.7060,0.6970,0.6910, & + & 0.6860,0.6830,0.6810,0.6800,0.6770,0.6730,0.6650,0.6560,0.6440, & + & 0.6320,0.6190,0.6080,0.5980,0.5900,0.5840,0.5790,0.5750,0.5710, & + & 0.5670,0.5600,0.5520,0.5420,0.5320,0.5210,0.5110,0.5020,0.4960, & + & 0.4910,0.4870,0.4840,0.4820,0.4790,0.4760,0.4730,0.4710,0.4700, & + & 0.4710,0.4740,0.4790,0.4860,0.4960,0.5080,0.5230,0.5400,0.5590, & + & 0.5790,0.5990,0.6200,0.6410,0.6630,0.6860,0.7110,0.7380,0.7680, & + & 0.7990,0.8330,0.8670,0.9010,0.9350,0.9680,1.0010,1.0330,1.0650, & + & 1.0970,1.1290,1.1620,1.1940,1.2270,1.2590,1.2910,1.3230,1.3550, & + & 1.3860,1.4170,1.4450,1.4730,1.4980,1.5220,1.5440,1.5650,1.5850, & + & 1.6060,1.6270,1.6470,1.6680,1.6900,1.7120,1.7350,1.7580,1.7810, & + & 1.8040,1.8270,1.8480,1.8700,1.8900,1.9100,1.9290,1.9470,1.9640, & + & 1.9810,1.9980,2.0150,2.0310,2.0490,2.0670,2.0850,2.1030,2.1210, & + & 2.1390,2.1550,2.1710,2.1870,2.2030,2.2180,2.2330,2.2480,2.2620, & + & 2.2740,2.2850,2.2960,2.3040,2.3120,2.3180,2.3220,2.3260,2.3300, & + & 2.3340,2.3390,2.3460,2.3540,2.3630,2.3720,2.3800,2.3860,2.3890, & + & 2.3890,2.3860,2.3810,2.3750,2.3680,2.3620,2.3560,2.3510,2.3460, & + & 2.3420,2.3390,2.3360,2.3340,2.3330,2.3310,2.3290,2.3260,2.3220, & + & 2.3180,2.3130,2.3100,2.3090,2.3080,2.3090,2.3080,2.3040,2.2960, & + & 2.2810,2.2600 /) + + real :: fsnowalb(965) = (/98.5210,98.6350,98.6700,98.7120,98.7790,98.8680,98.8450,98.9540,98.9610, & + & 99.0500,99.1250,99.2070,99.2820,99.2250,99.1880,99.2370,99.0990,99.1030, & + & 99.1020,99.1160,98.9910,98.9690,98.9930,98.9200,98.8790,98.7640,98.6570, & + & 98.5730,98.3690,98.2560,98.1060,97.7940,97.7790,97.7460,97.6460,97.4890, & + & 96.7700,97.0330,96.7160,96.6300,96.2220,95.9270,95.4900,95.5320,95.1660, & + & 94.6110,94.6240,94.5270,94.3100,93.3850,91.5990,91.2340,90.6630,90.1130, & + & 88.4780,86.1840,83.8720,82.0690,81.4550,83.1620,82.1010,84.7800,84.4700, & + & 82.9290,80.1690,75.9630,65.0150,67.5290,65.9010,64.8110,64.6430,65.2040, & + & 63.2540,65.8620,65.4960,64.2190,60.9990,44.1790,20.6330,13.9170,9.8620, & + & 8.5270,9.1870,9.9920,11.9100,13.6980,15.3350,16.5310,20.2890,20.9180, & + & 24.5220,28.3510,28.3660,29.4560,32.7410,35.2040,36.9860,38.7430,40.9240, & + & 37.5480,24.9240,15.3950,9.0950,6.2970,4.4940,3.6270,3.5060,3.6890, & + & 4.2150,5.0350,6.1990,8.1690,11.1840,13.9240,17.5410,20.1650,22.8010, & + & 26.3190,28.3300,27.8650,26.8480,25.0320,22.5830,19.9560,17.5520,15.3820, & + & 13.3700,11.7490,10.4760,9.5200,8.8530,8.3410,8.0240,7.7990,7.6360, & + & 7.4400,7.2310,6.9970,6.7890,6.5870,6.3510,6.1360,5.9150,5.7040, & + & 5.5560,5.4300,5.3400,5.2720,5.2760,5.2940,5.3000,5.3570,5.3790, & + & 5.2740,5.2360,5.1410,4.9320,4.6360,4.4280,4.0220,3.5180,2.8720, & + & 2.3400,1.7760,1.3320,0.9540,0.6770,0.5140,0.3780,0.3210,0.2640, & + & 0.2400,0.2350,0.2140,0.1900,0.1930,0.1480,0.1660,0.1710,0.1600, & + & 0.2100,0.1690,0.2100,0.2460,0.3040,0.3370,0.4280,0.5060,0.5860, & + & 0.6630,0.7190,0.7860,0.7970,0.8760,0.9160,1.0290,1.1780,1.2890, & + & 1.4920,1.7850,2.0290,2.3160,2.5890,2.7750,2.7460,2.5550,2.3160, & + & 2.1490,2.0570,2.0130,2.0030,2.0260,2.0420,2.0510,1.9820,1.9340, & + & 1.8640,1.7710,1.7000,1.6370,1.5850,1.5150,1.4780,1.4190,1.3810, & + & 1.3340,1.3090,1.3030,1.2660,1.2370,1.2160,1.2230,1.2200,1.1860, & + & 1.1730,1.1720,1.1580,1.1530,1.1390,1.1320,1.1310,1.1460,1.1470, & + & 1.1490,1.1580,1.1700,1.1670,1.1720,1.1860,1.2030,1.2210,1.2150, & + & 1.2240,1.2520,1.2750,1.2840,1.2860,1.2800,1.2810,1.2800,1.2860, & + & 1.2830,1.2730,1.2630,1.2360,1.2300,1.2140,1.1780,1.1460,1.1430, & + & 1.1320,1.0970,1.0680,1.0480,1.0360,1.0100,0.9780,0.9690,0.9610, & + & 0.9390,0.9150,0.9050,0.8920,0.8870,0.8840,0.8730,0.8470,0.8370, & + & 0.8390,0.8280,0.8120,0.7950,0.7960,0.7960,0.7840,0.7710,0.7670, & + & 0.7590,0.7550,0.7580,0.7540,0.7450,0.7320,0.7250,0.7280,0.7230, & + & 0.7210,0.7200,0.7140,0.7200,0.7180,0.7050,0.7050,0.7060,0.6910, & + & 0.6810,0.6860,0.6870,0.6690,0.6610,0.6700,0.6720,0.6560,0.6470, & + & 0.6550,0.6580,0.6510,0.6460,0.6410,0.6400,0.6400,0.6390,0.6380, & + & 0.6380,0.6320,0.6290,0.6350,0.6320,0.6360,0.6470,0.6460,0.6390, & + & 0.6350,0.6380,0.6430,0.6400,0.6280,0.6310,0.6540,0.6740,0.6740, & + & 0.6630,0.6560,0.6540,0.6560,0.6660,0.6730,0.6680,0.6630,0.6640, & + & 0.6700,0.6720,0.6710,0.6760,0.6880,0.6940,0.6880,0.6790,0.6820, & + & 0.6870,0.6880,0.6900,0.6930,0.6950,0.6990,0.7030,0.7020,0.7000, & + & 0.7000,0.7040,0.7110,0.7190,0.7290,0.7330,0.7330,0.7370,0.7390, & + & 0.7330,0.7290,0.7350,0.7440,0.7440,0.7290,0.7220,0.7310,0.7400, & + & 0.7380,0.7290,0.7270,0.7340,0.7420,0.7420,0.7350,0.7280,0.7210, & + & 0.7190,0.7200,0.7190,0.7110,0.7010,0.6960,0.6970,0.6950,0.6890, & + & 0.6820,0.6770,0.6730,0.6660,0.6560,0.6470,0.6430,0.6460,0.6500, & + & 0.6490,0.6420,0.6310,0.6250,0.6230,0.6200,0.6140,0.6090,0.6090, & + & 0.6130,0.6090,0.5980,0.5860,0.5770,0.5800,0.5820,0.5780,0.5660, & + & 0.5590,0.5610,0.5680,0.5730,0.5730,0.5700,0.5600,0.5490,0.5430, & + & 0.5410,0.5400,0.5380,0.5360,0.5310,0.5270,0.5190,0.5090,0.5050, & + & 0.5060,0.5080,0.5080,0.5030,0.4990,0.4970,0.4980,0.5000,0.4980, & + & 0.4910,0.4800,0.4720,0.4750,0.4840,0.4910,0.4890,0.4770,0.4640, & + & 0.4640,0.4720,0.4780,0.4740,0.4630,0.4510,0.4460,0.4520,0.4590, & + & 0.4620,0.4600,0.4580,0.4540,0.4520,0.4520,0.4560,0.4610,0.4640, & + & 0.4610,0.4530,0.4450,0.4350,0.4330,0.4430,0.4580,0.4660,0.4630, & + & 0.4560,0.4470,0.4460,0.4500,0.4530,0.4510,0.4420,0.4330,0.4320, & + & 0.4350,0.4450,0.4530,0.4630,0.4700,0.4730,0.4740,0.4750,0.4770, & + & 0.4770,0.4770,0.4750,0.4720,0.4700,0.4700,0.4710,0.4720,0.4700, & + & 0.4660,0.4580,0.4500,0.4460,0.4440,0.4490,0.4580,0.4660,0.4710, & + & 0.4740,0.4730,0.4720,0.4730,0.4760,0.4780,0.4790,0.4800,0.4810, & + & 0.4850,0.4880,0.4900,0.4920,0.4920,0.4910,0.4860,0.4830,0.4830, & + & 0.4860,0.4970,0.5090,0.5140,0.5160,0.5140,0.5110,0.5100,0.5120, & + & 0.5160,0.5180,0.5170,0.5120,0.5050,0.4990,0.4950,0.4980,0.5020, & + & 0.5040,0.5030,0.4980,0.4930,0.4900,0.4910,0.4980,0.5080,0.5130, & + & 0.5160,0.5180,0.5190,0.5190,0.5190,0.5200,0.5210,0.5240,0.5270, & + & 0.5300,0.5310,0.5260,0.5150,0.5060,0.4990,0.4950,0.4960,0.5000, & + & 0.5060,0.5120,0.5150,0.5160,0.5140,0.5110,0.5090,0.5070,0.5070, & + & 0.5090,0.5110,0.5110,0.5080,0.5020,0.4930,0.4860,0.4820,0.4820, & + & 0.4860,0.4920,0.4960,0.4960,0.4930,0.4880,0.4820,0.4790,0.4760, & + & 0.4760,0.4760,0.4760,0.4760,0.4780,0.4800,0.4810,0.4820,0.4830, & + & 0.4820,0.4800,0.4770,0.4740,0.4700,0.4670,0.4650,0.4630,0.4610, & + & 0.4580,0.4550,0.4540,0.4550,0.4600,0.4680,0.4770,0.4850,0.4910, & + & 0.4920,0.4890,0.4830,0.4750,0.4690,0.4650,0.4650,0.4670,0.4710, & + & 0.4760,0.4800,0.4820,0.4820,0.4810,0.4790,0.4760,0.4750,0.4750, & + & 0.4750,0.4750,0.4750,0.4740,0.4720,0.4710,0.4700,0.4700,0.4720, & + & 0.4720,0.4720,0.4700,0.4660,0.4630,0.4590,0.4570,0.4560,0.4570, & + & 0.4600,0.4640,0.4700,0.4760,0.4800,0.4830,0.4820,0.4810,0.4790, & + & 0.4770,0.4750,0.4750,0.4740,0.4740,0.4750,0.4770,0.4790,0.4790, & + & 0.4760,0.4710,0.4640,0.4590,0.4560,0.4560,0.4570,0.4570,0.4550, & + & 0.4520,0.4500,0.4510,0.4550,0.4610,0.4660,0.4700,0.4700,0.4700, & + & 0.4680,0.4650,0.4620,0.4580,0.4540,0.4510,0.4490,0.4490,0.4510, & + & 0.4550,0.4620,0.4710,0.4800,0.4860,0.4870,0.4820,0.4710,0.4590, & + & 0.4490,0.4420,0.4410,0.4430,0.4450,0.4450,0.4430,0.4410,0.4420, & + & 0.4460,0.4510,0.4540,0.4560,0.4530,0.4480,0.4440,0.4420,0.4450, & + & 0.4510,0.4590,0.4670,0.4750,0.4800,0.4820,0.4830,0.4830,0.4790, & + & 0.4720,0.4650,0.4590,0.4560,0.4570,0.4620,0.4680,0.4750,0.4810, & + & 0.4860,0.4910,0.4950,0.4990,0.5000,0.4970,0.4920,0.4850,0.4790, & + & 0.4730,0.4670,0.4570,0.4440,0.4280,0.4120,0.3980,0.3890,0.3850, & + & 0.3830,0.3800,0.3710,0.3580,0.3410,0.3250,0.3130,0.3070,0.3060, & + & 0.3090,0.3090,0.3050,0.2980,0.2890,0.2790,0.2680,0.2600,0.2530, & + & 0.2470,0.2430,0.2410,0.2420,0.2430,0.2460,0.2490,0.2510,0.2500, & + & 0.2460,0.2400,0.2340,0.2280,0.2250,0.2220,0.2190,0.2160,0.2140, & + & 0.2140,0.2170,0.2230,0.2310,0.2400,0.2480,0.2560,0.2650,0.2760, & + & 0.2900,0.3090,0.3300,0.3510,0.3710,0.3860,0.3970,0.4050,0.4120, & + & 0.4220,0.4330,0.4470,0.4640,0.4820,0.5010,0.5220,0.5450,0.5690, & + & 0.5930,0.6160,0.6350,0.6500,0.6600,0.6670,0.6710,0.6740,0.6790, & + & 0.6890,0.7050,0.7260,0.7510,0.7770,0.7980,0.8140,0.8230,0.8270, & + & 0.8310,0.8380,0.8490,0.8660,0.8840,0.9010,0.9150,0.9260,0.9340, & + & 0.9420,0.9520,0.9640,0.9790,0.9940,1.0090,1.0210,1.0320,1.0390, & + & 1.0440,1.0480,1.0520,1.0590,1.0690,1.0830,1.0990,1.1150,1.1260, & + & 1.1320,1.1330,1.1330,1.1350,1.1430,1.1560,1.1700,1.1790,1.1820, & + & 1.1780,1.1690,1.1570,1.1490,1.1460,1.1500,1.1610,1.1750,1.1880, & + & 1.1960,1.1940,1.1820,1.1630,1.1420,1.1230,1.1120,1.1100,1.1170, & + & 1.1280,1.1350,1.1360,1.1270,1.1130,1.0990,1.0920,1.0950,1.1060, & + & 1.1200,1.1280,1.1280,1.1170,1.0990,1.0770,1.0580,1.0440,1.0400, & + & 1.0440,1.0500,1.0660,1.0760,1.0790,1.0700,1.0520,1.0340,1.0200, & + & 1.0100,1.0100 /) + + + if (mypid == NROOT) then + + !snowalbedos(:) = 0.25*(fsnowalb(:)+2.0*msnowalb(:)+csnowalb(:)) !assume mostly med-grain + + wmin = const/(starbbtemp*36.841361) !Wavelength where exponential term is <=1.0e-16 + lwmin = log10(wmin) + + hinge = log10(7.5e-7) !We care about amounts above and below 0.75 microns + dl1 = (hinge-lwmin)/1024.0 + dl2 = (-4-hinge)/1024.0 + + do k=1,1024 + wv1(k) = 10**(lwmin+(k-1)*dl1) + wv2(k) = 10**(hinge+(k-1)*dl2) + enddo + do k=1,1024 + wvm1(k) = (1.0e6 * wv1(k))**5 + wvm2(k) = (1.0e6 * wv2(k))**5 + enddo + +! const1 = 2*planckh*(cc**2) + const2 = const/starbbtemp + + do k=1,1024 !Compute the Planck function + bb1(k) = 1.0/wvm1(k) * 1.0/(exp(const2/wv1(k))-1) !const1/wv1(k)**5 + bb2(k) = 1.0/wvm2(k) * 1.0/(exp(const2/wv2(k))-1) +! write(nud,*) wv1(k),bb1(k),wv2(k),bb2(k) + enddo !The scaling and units don't actually matter, because we're going to normalize + + do k=1,965 !Compute the Planck function for the wavelengths at which we have albedo data + bb3(k) = 1.0/(wavelengths(k))**5 * 1.0/(exp(1.0e6*const2/wavelengths(k))-1) + enddo + a1 = 0.0 + a2 = 0.0 + do k=1,41 !Compute insolation-weighted albedo below 0.75 microns + a1 = a1 + 0.5*(bb3(k)*fsnowalb(k)+bb3(k+1)*fsnowalb(k+1))* & + & 1.0e-6*(wavelengths(k+1)-wavelengths(k)) + enddo + do k=42,964 !Compute insolation-weighted albedo above 0.75 microns + a2 = a2 + 0.5*(bb3(k)*fsnowalb(k)+bb3(k+1)*fsnowalb(k+1)))* & + & 1.0e-6*(wavelengths(k+1)-wavelengths(k)) + enddo + + z1 = 0.0 + z2 = 0.0 + do k=1,1023 !Do a trapezoidal integration above and below 0.75 microns + z1 = z1 + 0.5*(bb1(k)+bb1(k+1))*(wv1(k+1)-wv1(k)) + z2 = z2 + 0.5*(bb2(k)+bb2(k+1))*(wv2(k+1)-wv2(k)) + enddo + z1 = z1 + 0.5*(bb1(1024)+bb2(1))*(wv2(1)-wv1(1024)) + + a1 = 0.01*a1/z1 !Percent -> Decimal; normalization + a2 = 0.01*a2/z2 + + znet = z1+z2 + + z1 = z1/znet + z2 = 1.0-z1 + + zsolar1 = z1 + zsolar2 = z2 + + write(nud,*) "Energy fraction below 0.75 microns:",zsolar1 + write(nud,*) "Energy fraction above 0.75 microns:",zsolar2 + + zsolars(1) = zsolar1 + zsolars(2) = zsolar2 + + dsnowalb(1) = a1 + dsnowalb(2) = a2 + + write(nud,*) "Snow albedo below 0.75 microns:",dsnowalb(1) + write(nud,*) "Snow albedo above 0.75 microns:",dsnowalb(2) + write(nud,*) "Overall snow albedo:",z1*dsnowalb(1)+z2*dsnowalb(2) + +! call put_restart_array("zsolars",zsolars,2,2,1) + + endif + + + + call mpbcrn(zsolars,2) + call mpbcr(zsolar1) + call mpbcr(zsolar2) + call mpbcrn(dsnowalb,2) + + call mpputgp('zsolars',zsolars,2,1) + call mpputgp('dsnowalb',dsnowalb,2,1) + + end subroutine solarini + ! ================= ! SUBROUTINE RADINI ! ================= @@ -145,7 +719,7 @@ & ,iyrbp,nswr,nlwr,nfixed,fixedlon,slowdown & & ,a0o3,a1o3,aco3,bo3,co3,toffo3,o3scale,newrsc & & ,nsol,nclouds,nswrcl,nrscat,rcl1,rcl2,acl2,clgray,tpofmt & - & ,acllwr,tswr1,tswr2,tswr3,th2oc,dawn + & ,acllwr,tswr1,tswr2,tswr3,th2oc,dawn,starbbtemp,nstartemp ! ! namelist parameter: ! @@ -306,6 +880,37 @@ call mpbci(nswrcl) call mpbci(nclouds) + call mpbcr(starbbtemp) + call mpbci(nstartemp) + +! +! determine stellar parameters +! + if (nstartemp > 0) then + if (nrestart > 0.) then + call mpgetgp('zsolars',zsolars,2,1) + call mpgetgp('dsnowalb',dsnowalb,2,1) + if (mypid == NROOT) then +! call get_restart_array("zsolars",zsolars,2,2,1) + zsolar1 = zsolars(1) + zsolar2 = zsolars(2) + write(nud,*) "Read zsolar1 from restart: ",zsolar1 + write(nud,*) "Read zsolar2 from restart: ",zsolar2 + write(nud,*) "Read snow albedo <0.75 um from restart: ",dsnowalb(1) + write(nud,*) "Read snow albedo >0.75 um from restart: ",dsnowalb(2) + endif + call mpbcr(zsolar1) + call mpbcr(zsolar2) + call mpbcrn(dsnowalb) +! call mpputgp( + else + call solarini + endif + else + call mpbcr(zsolar1) + call mpbcr(zsolar2) + call mpbcrn(dsnowalb) + endif ! ! determine orbital parameters ! @@ -395,7 +1000,8 @@ real, allocatable :: zprf11(:,:) real, allocatable :: zprf12(:,:) real, allocatable :: zcc(:,:) - real, allocatable :: zalb(:) + real, allocatable :: zalb1(:) + real, allocatable :: zalb2(:) real, allocatable :: zdtdte(:,:) ! ! cpu time estimates @@ -429,9 +1035,11 @@ if(ndiagcf > 0) then allocate(zcc(NHOR,NLEP)) - allocate(zalb(NHOR)) + allocate(zalb1(NHOR)) + allocate(zalb2(NHOR)) zcc(:,:)=dcc(:,:) - zalb(:)=dalb(:) + zalb1(:) = dalb1(:) + zalb2(:) = dalb2(:) dcc(:,:)=0. if(nswr==1) call swr dclforc(:,1)=dswfl(:,NLEP) @@ -439,8 +1047,10 @@ dclforc(:,5)=dfu(:,1) dclforc(:,6)=dfu(:,NLEP) dcc(:,:)=zcc(:,:) - dalb(:)=zalb(:) - deallocate(zalb) + dalb1(:) = zalb1(:) + dalb2(:) = zalb2(:) + deallocate(zalb1) + deallocate(zalb2) end if ! @@ -712,6 +1322,8 @@ ! ! no PUMA variables are used ! + call mpputgp('zsolars',zsolars,2,1) + if(mypid == NROOT .and. ntime == 1) then write(nud,*)'******************************************' write(nud,*)' CPU usage in RADSTEP (ROOT process only): ' @@ -941,8 +1553,6 @@ ! 0) define local parameters and arrays ! parameter(zero=1.E-6) ! if insolation < zero : fluxes=0. - parameter(zsolar1=0.517) ! spectral partitioning 1 (wl < 0.75mue) - parameter(zsolar2=0.483) ! spectral partitioning 2 (wl > 0.75mue) parameter(zbetta=1.66) ! magnification factor water vapour parameter(zmbar=1.9) ! magnification factor ozon parameter(zro3=2.14) ! ozon density (kg/m**3 STP) @@ -1322,23 +1932,26 @@ ! make upward R ! -! -! This would be the place to add spectral albedo dependence for snow -! and ice. Or at least one of the options. The other option is to -! keep track of two dalb arrays, and when dalb is computed, also -! compute a NIR version. The latter is more computationally expensive, -! so we will check first if we can use a simple multiplier. -! -! - zra1s(:)=dalb(:) - zra2s(:)=dalb(:) +! Currently: we use the same albedo for both spectral ranges. + + + zra1s(:)=dalb1(:) + zra2s(:)=dalb2(:) ! ! set albedo for the direct beam (for ocean use ECHAM3 param) ! - dalb(:)=dls(:)*dalb(:)+(1.-dls(:))*dicec(:)*dalb(:) & + dalb1(:)=dls(:)*dalb1(:)+(1.-dls(:))*dicec(:)*dalb1(:) & & +(1.-dls(:))*(1.-dicec(:))*AMIN1(0.05/(zmu0(:)+0.15),0.15) - zra1(:)=dalb(:) - zra2(:)=dalb(:) + dalb2(:)=dls(:)*dalb2(:)+(1.-dls(:))*dicec(:)*dalb2(:) & + & +(1.-dls(:))*(1.-dicec(:))*AMIN1(0.05/(zmu0(:)+0.15),0.15) + + dalb(:) = zsolars(1)*dalb1(:) + zsolars(2)*dalb2(:) + + zra1(:)=dalb1(:) + zra2(:)=dalb2(:) + +! Ice-free ocean albedo is min(0.05/(phi+0.15), 0.15)--reflection and scattering is higher at low phi + endwhere do jlev=NLEV,1,-1 where(losun(:))