C+ C NAME: C Dust C PURPOSE: C Test program for the subroutine DustAsymmetry which transforms the Helios C photometer brightnesses into the plane of symmetry of the zodiacal dust cloud. C CATEGORY: C Test: Helios C CALLING SEQUENCE: program Dust C INPUTS: C Unnormalized Helios photometer, containing a single color and filter. C OUTPUTS: C To screen C CALLS: C BadR4, pInfR4, iHOSInfo, iHOSRead, iSetFileSpec, iGetFileSpec, ArrR4Copy C DustAsymmetry, AskYN, iSetTermAndPlot, iSwitchGraphicsOn, iSwitchGraphicsOff C ClearWindow, USRDF, PlotCurve, Hardcopy, AskI4, AskLogG C INCLUDE: include 'filparts.h' include 'hos_e9.h' C PROCEDURE: C Only filters 1 through 4 are analyzed. Filter 5 (pB) is ignored. C MODIFICATION HISTORY: C MAY-2000, Paul Hick (UCSD; pphick@ucsd.edu); documentation C- parameter (MXC=2500) parameter (NS =8) ! Start sector (modified) parameter (NST=1-NS) ! End sector (modified) parameter (NSTOT=NS-NST+1) ! Total # sectors integer ICNT(2) /2*0/ integer NN(MXC) integer PP(MXC) integer CC(MXC) integer FF(MXC) integer Lst(3,11) real TT(MXC) real RR(MXC) real LL(MXC) real ZZraw(NST:NS,MXC) real ZZsym(NST:NS,MXC) real ZZbck(NST:NS,MXC) real T1(MXC) real Z1(MXC) real T2(MXC) real Z2(MXC) character cFileIn*80 character cDir*40 character cSay*4 /'Dust'/ logical bSym1 logical bSym2 logical bTmp logical bOneCC logical bOneFF bOneCC(kC) = kC .eq. HOS__C_1 .or. kC .eq. HOS__C_2 .or. kC .eq. HOS__C_3 bOneFF(kF) = kF .eq. HOS__F_1 .or. kF .eq. HOS__F_2 .or. kF .eq. HOS__F_3 .or. kF .eq. HOS__F_4 .or. kF .eq. HOS__F_5 Bad = BadR4() Big = pInfR4() if (iHOSInfo(0,cFileIn,iRecl,iSc,iYr,TMin,TMax,kPCFall,Lst) .eq. 0) & call Say(cSay,'W','StopPlay','Empty input file or read error: '//cFileIn) if (iand(kPCFall,HOS__NONORM) .eq. 0) & call Say(cSay,'W','StopPlay','Data are already normalized: '//cFileIn) kC = iand(kPCFall,HOS__C_ALL) if (.not. bOneCC(kC)) call Say(cSay,'W','StopPlay','Only one color allowed: '//cFileIn) !------- ! Read 15/30 degree data records : ICNT(1) = # 16/31 deg records ! We are not interested in pB data. kPCF = kPCFAll kPCF = kPCF-iand(kPCF,HOS__P_ALL)+HOS__P_12+HOS__SWAP_SECT ! Phot 1&2; swap sectors kPCF = kPCF-iand(kPCF,HOS__F_ALL)+HOS__F_POL+HOS__CLR ! Exlude pB call HOSRead(kPCF,cFileIn,MXC,ICNT(1),NN,TT,PP,CC,FF,RR,LL,NSTOT,ZZraw) kF = iand(kPCF,HOS__F_ALL) if (.not. bOneFF(kF)) call Say(cSay,'W','StopPlay','Only one filter allowed: '//cFileIn) if (kF .eq. HOS__F_4) then !------- ! Read 90 deg data records : ICNT(2)-ICNT(1) = # 90 deg records kPCF = kPCFAll kPCF = kPCF-iand(kPCF,HOS__P_ALL)+HOS__P_3 ! Phot. 3 I = ICNT(1)+1 call HOSRead(kPCF,cFileIn,MXC,ICNT(2),NN(I),TT(I),PP(I),CC(I),FF(I),RR(I),LL(I),NSTOT,ZZraw(NST,I)) end if ICNT(2) = ICNT(1)+ICNT(2) print *, 'S/C = ', iSc, ' Year = ', iYr !------- ! Split off directory name from input file name I = iSetFileSpec(cFileIn) I = iGetFileSpec(0,FIL__DIRECTORY,cDir) I = iGetFileSpec(FIL__NAME,0,cFileIn) ! ========== Transform data to plane of symmetry of zodiacal dust cloud call ArrR4Copy(NSTOT*ICNT(2),ZZraw,ZZsym) call DustAsymmetry(2,iSc,ICNT(2),PP,LL,NSTOT,ZZsym) call AskYN('Backward transformation$no',bTmp) call AskI4('Photometer: 16, 31',IP) IS1 = 1 call AskI4('Select 1st sector$1$32',IS1) if (IS1 .gt. 16) IS1 = IS1-32 if (bTmp) then call ArrR4Copy(NSTOT*ICNT(2),ZZsym,ZZbck) call DustAsymmetry(3,iSc,ICNT(2),PP,LL,NSTOT,ZZbck) !------- ! ZZraw is the input intensity; ZZsym is transformed to the plane ! of symmetry; ZZbck is obtained by transforming ZZsym back and ! should be the same as ZZraw do I=1,ICNT(1) if (PP(I) .eq. IP .and. ZZraw(IS1,I) .ne. Bad) & print *, ZZraw(IS1,I), ZZsym(IS1,I)-ZZraw(IS1,I), ZZbck(IS1,I)-ZZraw(IS1,I) end do call Say('Dust','I','StopPlay','Done') end if call AskYN('Plot corrected data for 1st sector$no',bSym1) IS2 = 1 call AskI4('Select 2nd sector$1$32',IS2) if (IS2 .gt. 16) IS2 = IS2-32 call AskYN('Plot corrected data for 2nd sector$no',bSym2) if (.not. bSym1 .and. .not. bSym2) then do I=1,ICNT(1) if (PP(I) .eq. IP .and. ZZraw(IS1,I) .ne. Bad .and. ZZraw(IS2,I) .ne. Bad) & print *, ZZraw(IS1,I),ZZraw(IS2,I),(ZZraw(IS2,I)-ZZraw(IS1,I))/ZZraw(IS2,I) end do else if (bSym1 .and. .not. bSym2) then do I=1,ICNT(1) if (PP(I) .eq. IP .and. ZZsym(IS1,I) .ne. Bad .and. ZZraw(IS2,I) .ne. Bad) & print *, ZZsym(IS1,I),ZZraw(IS2,I),(ZZraw(IS2,I)-ZZsym(IS1,I))/ZZraw(IS2,I) end do else if (.not. bSym1 .and. bSym2) then do I=1,ICNT(1) if (PP(I) .eq. IP .and. ZZraw(IS1,I) .ne. Bad .and. ZZsym(IS2,I) .ne. Bad) & print *, ZZraw(IS1,I),ZZsym(IS2,I),(ZZsym(IS2,I)-ZZraw(IS1,I))/ZZsym(IS2,I) end do else if (bSym1 .and. bSym2) then do I=1,ICNT(1) if (PP(I) .eq. IP .and. ZZsym(IS1,I) .ne. Bad .and. ZZsym(IS2,I) .ne. Bad) & print *, ZZsym(IS1,I),ZZsym(IS2,I),(ZZsym(IS2,I)-ZZsym(IS1,I))/ZZsym(IS2,I) end do end if ZMin = Big ZMax = -Big N1 = 0 do I=1,ICNT(1) if (PP(I) .eq. IP .and. ZZsym(IS1,I) .ne. Bad) then N1 = N1+1 T1(N1) = TT(I) if ( bSym1) Z1(N1) = ZZsym(IS1,I) if (.not. bSym1) Z1(N1) = ZZraw(IS1,I) ZMin = min(ZMin,Z1(N1)) ZMax = max(ZMax,Z1(N1)) end if end do N2 = 0 do I=1,ICNT(1) if (PP(I) .eq. IP .and. ZZsym(IS2,I) .ne. Bad) then N2 = N2+1 T2(N2) = TT(I) if ( bSym2) Z2(N2) = ZZsym(IS2,I) if (.not. bSym2) Z2(N2) = ZZraw(IS2,I) ZMin = min(ZMin,Z2(N2)) ZMax = max(ZMax,Z2(N2)) end if end do if (ZMin .eq. Big) call Say('Dust','E','Stop','No data in file') ! ========== Display results call AskYN('Graphics output$no',bTmp) if (bTmp) then I = iSetTermAndPlot(-1,1) I = iSwitchGraphicsOn() call ClearWindow call USRDF(100.,100.,900.,700., TMin,ZMin,TMax,ZMax) call PlotCurve(1,1.,N1,T1,Z1) call PlotCurve(1,1.,N2,T2,Z2) call USRDF(0.,0.,1000.,700., 0.,0.,100.,100.) call AskLogG(30.,1.,'Hardcopy$no',bTmp) if (bTmp) call Hardcopy(2) I = iSwitchGraphicsOff() end if end