subroutine sums(slon,slat,mlon,mlat,ae,etos) 


	real mlat,mlon 
	common/grid/gridlon(128),gridlat(64),clat(64),mask(128,64) 
     1		      ,re(3,128,64),dlat(64) 
	common/clouds/ic(128,64) 
        common/contrib/contrib(128,64)
	common/ice/ice(128,64) 
	dimension rs(3),rm(3) 
	double precision dsind, dcosd, dasind, dacosd, datan2d
	external sind, cosd, dsind, dcosd, dasind, dacosd
        external datan2d, acosd, asind, atan2d
	



	rs(3)=sind(slat) 
	rm(3)=sind(mlat) 
	rs(1)=cosd(slat)*cosd(slon) 
	rm(1)=cosd(mlat)*cosd(mlon) 
	rs(2)=cosd(slat)*sind(slon) 
	rm(2)=cosd(mlat)*sind(mlon) 
	dotsm=dot(rs,rm) 
	ae=0. 
	etos=0.



	do j=1,64 
	   wt=clat(j)*dlat(j) 
	   do i=1,128 
	      dotes=dot(re(1,i,j),rs)
	      if (dotes.ge.0) then
	         if (dotes.ge.0.999999) dotes=0.999999 
	         zenith=acosd(dotes) 
	         itype1=mask(i,j)
	         icl=ic(i,j)
	         if (icl.ge.95) then
	            itype=12 
	         else if (icl.ge.50) then 
		       if (itype1.ge.2.and.itype1.le.4) then 
		          itype=10 
	             else if(itype1.eq.1) then 
		          itype=9 
		       else 
		          itype=11 
		       endif 
	         else if (icl.ge.5) then 
		       if (itype1.ge.2.and.itype1.le.4) then 
		          itype=7 
		       else if(itype1.eq.1) then 
		          itype=6 
		       else 
		          itype=8 
		       endif 
	         else 
	             itype=itype1 
	             if (ice(i,j).ne.0) itype=3 
	         endif
 
c******************************* 
c                 write(*,*) 'Vixa 2'
c******************************* itype=2 

	         asav=albedo(zenith,itype)
	         ae=ae+asav*dotes*wt 
	         dotem=dot(re(1,i,j),rm) 

	         if (dotem.ge.0.) then 
		         if (dotem.ge.0.999999) dotem=0.999999 
		         contrib(i,j)=1
	               view=acosd(dotem) 
	               azimuth=dotem*dotes-dotsm 
	               azimuth=azimuth/sqrt(1.-dotem**2) 
	               azimuth=azimuth/sqrt(1.-dotes**2) 
	               if (abs(azimuth).le.0.999999) then 
		               azimuth=acosd(azimuth) 
	               else if (azimuth.ge.0.999999) then
	                     azimuth=0. 
	               else
	                     azimuth=180. 
	               endif 
	               asav=asav*anifac(azimuth,view,zenith,itype) 
	               asym=albedo(view,itype)*anifac(azimuth,zenith,view,itype) 
	               asym=(asym+asav)/2. 
	               etos=etos+asym*wt*dotem*dotes 
	         endif 
	      endif 
	   enddo 
	enddo 

	ae=ae/64. 
	etos=etos/64. 


	return 
	end