c CODE_TO_DONES/mcnewd_oppo.for c**************************************************************** program mcnewd_oppo c**************************************************************** c driver for calculating of photometric properties c uses mcnewd_sub.for for the calculation of scattering c from a single position-file / sun direction) c several observing directions corresponding to c a fixed observing latitude and different phase-angles c output of results via wwwfile created in mcnewd_sub.for c use mc_wread_f.pro for reading the output c**************************************************************** c HS 7/12/99 - 04/01/02 implicit double precision (a-h,o-z) include 'mcnew_main.inc' parameter(nnt=1000) character*60 posfile,parfile,outfile,epsfile,wwwfile double precision testfii(nnt),testlat(nnt) double precision testw(nnt),testd(nnt),testn(nnt) double precision testm(nnt,10) double precision lonsun,mu0 double precision iscat_mean,npot_mean,nbox_mean character*8 date character*10 time character*5 zone integer*4 values(8) character*12 smet double precision ar,ra,apu1,apu2,apu3,apu4 ar=0.01745329251994329547d0 ra=1.d0/ar c posfile = posdat-file containing particle positions c parfile = parameter-file (junk -> defaults) c outfile = output of weights + directions (junk) c epsfile = output of weights + fiilat (junk) c wwwfile = output of weights c sredu = reduction factor for particle radius (sredu < 1) c zredu = reduction factor for vertical coordinates (zredu > 1) c idum = seed c******************************************************** c READ INPUT VALUES c******************************************************** c particle positions from POSFILE c BSUN,LONSUN = elevation and longitude of sun c TESTLAT,TESTFII = table of observer positions c OUTPUT: c testw,testn,testd,testm c to convert testw to I/F: c multiply by sin(bsun)/nphot/sin(bobs) read(*,*) posfile read(*,*) bsun,lonsun read(*,*) ntest c testfii intepreted as alpha c cos_alpha=cos(dlobs*ar)*cos(bsun*ar)*cos(bobs*ar)+ c sin(ar*bobs)*sin(bsun*ar) do 10 i=1,ntest read(*,*) testlat(i),testfii(i) apu=testfii(i) apu1=dcos(testfii(i)*ar)-dsin(ar*bsun)*dsin(ar*testlat(i)) apu2=dcos(ar*bsun)*dcos(ar*testlat(i)) apu3=apu1/apu2 apu4=apu3 apu3=dacos(apu3)*ra if(abs(testfii(i)).le.1d-6) apu3=0.d0 testfii(i)=lonsun+apu3 write(*,*) 'lat,fii,alpha',testlat(i),testfii(i),apu 10 continue read(*,*) nphot,mx,my read(*,*) parfile read(*,*) outfile read(*,*) epsfile read(*,*) wwwfile read(*,*) method smet='unknown' if(method.eq.1 ) smet='LAMBERT-element ' if(method.eq.11) smet='LAMBERT-sphere-center ' if(method.eq.12) smet='LAMBERT-sphere-point ' if(method.eq.2 ) smet='ISOTROPIC-sphere' if(method.eq.3 ) smet='HG-sphere-center' if(method.eq.4 ) smet='HG-sphere-point' if(method.eq.13) smet='HG2-sphere-center' if(method.eq.14) smet='HG2-sphere-point' if(method.eq.5 ) smet='EULER-sphere-center' if(method.eq.6 ) smet='EULER-sphere-point' if(method.eq.7 ) smet='CALLISTO-sphere-center' if(method.eq.8 ) smet='CALLISTO-sphere-point' if(method.eq.9 ) smet='POWER-sphere-center' if(method.eq.10) smet='POWER-sphere-point' read(*,*) ppar1,ppar2,ppar3 c icheck=-1 -> no check output read(*,*) icheck read(*,*) idum read(*,*) sredu,zredu read(*,*) sizemin read(*,*) s20,f20 c----------------------------------------- c Saturn-shine: c iglobe=0 -> SUN c iglobe=1 -> SATURN c iglobe=2 -> SATURN-SOUTH c iglobe=3 -> SATURN-NORTH c abox = calculation distance c rsaturn = Saturn radius c nlat_sat = number of latitude bins c nfii_sat = number of longitude bins c model_sat = saturn model c = 1 -> Lambert c = 2 -> RED (DONES ET ALL) c = 3 -> BLUE (DONES ET ALL) rsaturn=60.4 read(*,*) iglobe read(*,*) abox read(*,*) nlat_sat,nfii_sat read(*,*) model_sat if(model_sat.eq.1) smodel='Lambert-planet' if(model_sat.eq.2) smodel='RED 0.5*zones+0.5*belts' if(model_sat.eq.3) smodel='BLUE 0.5*zones+0.5*belts' if(model_sat.eq.4) smodel='RED+BLUE 0.5*zones+0.5*belts' if(model_sat.eq.5) smodel='ZONES 0.75*blue+0.25*red' if(model_sat.eq.6) smodel='BELTS 0.75*blue+0.25*red' write(*,*) '==========================================' write(*,*) 'CODE_TO_DONES/mcnewd_oppo.for (04.01.2002)' call date_and_time (date, time, zone, values) write(*,999) values(1),values(2),values(3),values(5),values(6) 999 format(1x,i4,'-',i2,'-',i2,2x,i2,':',i2,':',i2) write(*,*) '==========================================' write(*,*) 'POSFILE= ',posfile write(*,*) 'BSUN, LONSUN=',bsun,lonsun write(*,*) 'NTEST ',ntest write(*,*) 'NPHOT,MX,MY= ',nphot,mx,my write(*,*) 'METHOD = ',method,' ',smet write(*,*) 'PPAR 1,2,3 = ',ppar1,ppar2,ppar3 write(*,*) 'SEED = ',idum write(*,*) 'parfile,outfile,epsfile,wwwfile=' write(*,*) '-> ',parfile write(*,*) '-> ',outfile write(*,*) '-> ',epsfile write(*,*) '-> ',wwwfile write(*,*) 'sredu,zredu',sredu,zredu write(*,*) 'sizemin ',sizemin write(*,*) 's20,f20 ',s20,f20 if(iglobe.eq.0) write(*,*) 'ILLUMINATION BY SUN' if(iglobe.eq.1) write(*,*) 'ILLUMINATION BY SATURN' if(iglobe.eq.2) write(*,*) 'ILLUMINATION BY SATURN-SOUTH' if(iglobe.eq.3) write(*,*) 'ILLUMINATION BY SATURN-NORTH' if(iglobe.ge.1) then write(*,*) 'MODEL_SAT=',model_sat,smodel write(*,*) 'DISTANCE =',abox,abox/rsaturn write(*,*) 'NLAT,NFII=',nlat_sat,nfii_sat endif write(*,*) 'FIIOBS intepreted as alpha' write(*,*) '==========================================' sunfii=lonsun sunlat=bsun c Saturn-shine? nphot_tod=nphot flux_sat=0. if(iglobe.ge.1) then call saturn_illumination + (sunlat,sunfii,nphot_ori,nphot_tod,flux_sat) nphot=nphot_tod endif write(*,*) 'flux_sat,nphot',flux_sat,nphot c MC-calculation call tim(cpu1) call mcnew_sub(nphot,posfile,parfile,outfile, + sunfii,sunlat, + ntest,testfii,testlat,testw,testd,testn,testm, + idum,icheck,mx,my, + npart,npartori,npot_mean,nbox_mean,iscat_mean, + method,epsfile,wwwfile,ifor,sredu,tau,ndirect, + zredu,ppar1,ppar2,ppar3) call tim(cpu2) mu0=cos((90.-sunlat)*0.017453293) fdirect=1.*ndirect/nphot tau_path=-log(fdirect) tau_n=tau_path*mu0 cpu_comp=(cpu2-cpu1)/(1.*nphot/1000.) write(*,*) 'NORI,TAU = ',npartori,tau write(*,*) 'SREDU, ZREDU = ', sredu,zredu write(*,*) 'cpu, cpu_comp',cpu2-cpu1,cpu_comp write(*,*) 'npart,npot,nbox,iscat', + npart,npot_mean,nbox_mean,iscat_mean write(*,*) 'fdir,mu0.tau_p,tau_n', + fdirect,mu0,tau_path,tau_n write(*,*) 'sizemin ',sizemin write(*,*) 's20,f20 ',s20,f20 write(*,*) '==========================================' call date_and_time (date, time, zone, values) write(*,999) values(1),values(2),values(3),values(5),values(6) write(*,*) '==========================================' end