PROGRAM BROADENX c constants given D exponents c revised 31jan01 changed read statement for wledge c C TAPE5=INPUT C TAPE6=OUTPUT C TAPE21=SPECTRUM INPUT C TAPE22=SPECTRUM OUTPUT C THE MINIMUM DIMENSION OF H IS (NWL+19999+19999)*NMU C FOR FLUX SPECTRA NMU IS 1 C DIMENSION H(4000000) DIMENSION H(4100000) DIMENSION RED(20000),BLUE(20000) DIMENSION RED1(20000),BLUE1(20000),RED2(20000),BLUE2(20000) EQUIVALENCE (RED(1),RED1(1)),(BLUE(1),BLUE1(1)) CHARACTER*10 A,B CHARACTER*50 COMMENT REAL*8 LINDAT8(14) REAL*4 LINDAT4(28) DIMENSION XMU(20),QMU(40),WLEDGE(377),TITLE(74) REAL*8 TEFF,GLOG,TITLE,WBEGIN,RESOLU,XMU,WLEDGE,RATIO REAL*8 QMU CHARACTER*1 APLOT(101) DATA APLOT/101*' '/ LINOUT=300 REWIND 21 READ(21)TEFF,GLOG,TITLE,WBEGIN,RESOLU,NWL,IFSURF,NMU,XMU,NEDGE, 1(WLEDGE(iedge),iedge=1,nedge) WRITE(6,1010)TEFF,GLOG,TITLE 1010 FORMAT( 5H TEFF,F7.0,7H GRAV,F7.3/7H TITLE ,74A1) WRITE(6,1007)NMU,(XMU(IMU),IMU=1,NMU) 1007 FORMAT(I4,20F6.3) C IFSURF=3 FOR ROTATED SPECTRUM IF(IFSURF.EQ.3) NMU=1 RATIO=1.+1./RESOLU WEND=WBEGIN*RATIO**(NWL-1) WCEN=(WBEGIN+WEND)*.5 VSTEP=2.99792458D5/RESOLU WRITE(6,1005)WBEGIN,WEND,RESOLU,VSTEP,NWL 1005 FORMAT(2F14.5,F12.1,F12.5,I10) NMU1=NMU+1 NMU2=NMU+NMU C C SAMPLE CARDS RIGHT SHIFTED BY 1 C THESE ARE EVALUATED AT WLCEN AND HAVE CONSTANT C RESOLUTION OR RESOLVING POWER CGAUSSIAN 3.5 KM COMMENT FIELD CGAUSSIAN 100000. RESOLUTION COMMENT FIELD CGAUSSIAN 7. PM COMMENT FIELD CGAUSSIAN .01 CM-1 COMMENT FIELD CSINX/X 3.5 KM COMMENT FIELD CSINX/X 100000. RESOLUTION COMMENT FIELD CSINX/X 7. PM COMMENT FIELD CSINX/X .01 CM-1 COMMENT FIELD CRECT 7. PM COMMENT FIELD CRECT 3.5 KM COMMENT FIELD CRECT 100000. RESOLUTION COMMENT FIELD CREDT .01 CM-1 COMMENT FIELD CMACRO 2.0 KM COMMENT FIELD CPROFILE 5. POINTS COMMENT FIELD CRED .3 .1 .1 .1 .05 CBLUE .3 .1 .1 .1 .05 C POINTS ARE TABULATED AT THE SPACING OF THE COMPUTED SPECTRUM C THE CENTER IS THE FIRST POINT FOR EACH WING C TAKING THE CENTER ONLY ONCE, THE PROFILE SUMS TO 1. C C IN THESE PROFILES THE RESOLUTION OR RESOLVING POWER VARIES C THEY ARE LINEARLY INTERPOLATED FROM WLBEG TO WLEND C AT WLBEG AT WLEND CGAUSSIAN 3.5 KM 3. COMMENT FIELD CGAUSSIAN 100000. RESOLUTION 120000. COMMENT FIELD CGAUSSIAN 7. PM 7. COMMENT FIELD CGAUSSIAN .01 CM-1 .01 COMMENT FIELD CSINX/X 3.5 KM 3. COMMENT FIELD CSINX/X 100000. RESOLUTION 120000. COMMENT FIELD CSINX/X 7. PM 7. COMMENT FIELD CSINX/X .01 CM-1 .01 COMMENT FIELD CRECT 7. PM 7. COMMENT FIELD CRECT 3.5 KM 3. COMMENT FIELD CRECT 100000. RESOLUTION 120000. COMMENT FIELD CRECT .01 CM-1 .01 COMMENT FIELD CPROFILE 5. POINTS 5. COMMENT FIELD CRED .3 .1 .1 .1 .05 CBLUE .3 .1 .1 .1 .05 CRED .25 .125 .1 .1 .05 CBLUE .25 .125 .1 .1 .05 READ(5,1)A,X1,B,X2,COMMENT 1 FORMAT(A10,F10.2,A10,F10.2,A50) WRITE(6,2)A,X1,B,X2,COMMENT 2 FORMAT(1X,A10,F10.2,A10,F10.2,A50) FWHM=-1. IF(B.EQ.'PM ')FWHM=X1/WCEN/1000.*299792.458D0 IF(B.EQ.'KM ')FWHM=X1 IF(B.EQ.'RESOLUTION')FWHM=299792.458D0/X1 IF(B.EQ.'CM-1 ')FWHM=X1/(1.E7/WCEN)*299792.458D0 IF(FWHM.LT.0.)THEN WRITE(6,3) CALL EXIT ENDIF FWHM1=FWHM FWHM2=FWHM IF(X2.GT.0.)THEN IF(B.EQ.'PM ')FWHM2=X2/WCEN/1000.*299792.458D0 IF(B.EQ.'KM ')FWHM2=X2 IF(B.EQ.'RESO ')FWHM2=299792.458D0/X2 IF(B.EQ.'CM-1 ')FWHM2=X2/(1.E7/WCEN)*299792.458D0 ENDIF IF(A.EQ.'MACRO ')GO TO 10 IF(A.EQ.'GAUSSIAN ')GO TO 20 IF(A.EQ.'SINX/X ')GO TO 60 IF(A.EQ.'RECT ')GO TO 30 IF(A.EQ.'PROFILE ')GO TO 40 WRITE(6,3) 3 FORMAT(10H0BAD INPUT) CALL EXIT C C MACROTURBULENT VELOCITY IN KM 10 VMAC=X1 DO 11 I=1,20000 RED(I)=EXP(-(FLOAT(I-1)*VSTEP/VMAC)**2) IF(RED(I).LT.1.E-5)GO TO 12 11 CONTINUE 12 NPROF=I RED(1)=RED(1)/2. SUM=0. DO 13 I=1,NPROF 13 SUM=SUM+RED(I) SUM=SUM*2. DO 14 I=1,NPROF RED(I)=RED(I)/SUM 14 BLUE(I)=RED(I) GO TO 50 C C GAUSSIAN INSTRUMENTAL PROFILE HALF WIDTH IN KM FWHM 20 DO 21 I=1,20000 RED(I)=EXP(-(FLOAT(I-1)*VSTEP/FWHM*.8325546D0*2.)**2) IF(RED(I).LT.1.E-5)GO TO 22 21 CONTINUE 22 NPROF=I RED(1)=RED(1)/2. SUM=0. DO 23 I=1,NPROF 23 SUM=SUM+RED(I) SUM=SUM*2. DO 24 I=1,NPROF RED(I)=RED(I)/SUM 24 BLUE(I)=RED(I) IF(X2.EQ.0.)GO TO 50 DO 25 I=1,20000 RED2(I)=EXP(-(FLOAT(I-1)*VSTEP/FWHM2*.8325546D0*2.)**2) IF(RED2(I).LT.1.E-5)GO TO 26 25 CONTINUE 26 NPROF=I RED2(1)=RED2(1)/2. SUM=0. DO 27 I=1,NPROF 27 SUM=SUM+RED2(I) SUM=SUM*2. DO 28 I=1,NPROF RED2(I)=RED2(I)/SUM 28 BLUE2(I)=RED2(I) GO TO 50 C C SINX/X INSTRUMENTAL PROFILE HALF WIDTH IN KM FWHM C APODIZED BY EXP(-0.06*X**2) 60 RED(1)=0.5 DO 61 I=2,20000 X=(FLOAT(I-1)*VSTEP/FWHM*2.*1.8954942D0) RED(I)=SIN(X)/X*EXP(-0.06*X**2) IF(ABS(RED(I))+ABS(RED(I-1)).LT.1.E-5)GO TO 62 61 CONTINUE 62 NPROF=I SUM=0. DO 63 I=1,NPROF 63 SUM=SUM+RED(I) SUM=SUM*2. DO 64 I=1,NPROF RED(I)=RED(I)/SUM 64 BLUE(I)=RED(I) IF(X2.GT.0)GO TO 50 RED2(1)=0.5 DO 65 I=2,20000 X=(FLOAT(I-1)*VSTEP/FWHM2*2.*1.8954942D0) RED2(I)=SIN(X)/X*EXP(-0.06*X**2) IF(ABS(RED2(I))+ABS(RED2(I-1)).LT.1.D-5)GO TO 66 65 CONTINUE 66 NPROF=I SUM=0. DO 67 I=1,NPROF 67 SUM=SUM+RED2(I) SUM=SUM*2. DO 68 I=1,NPROF RED2(I)=RED2(I)/SUM 68 BLUE2(I)=RED2(I) GO TO 50 C C RECTANGULAR INSTRUMENTAL PROFILE HALF WIDTH IN KM FWHM 30 XRECT=FWHM/2./VSTEP NRECT=XRECT+1.5 NPROF=NRECT DO 31 I=1,NPROF 31 RED(I)=1. RED(NPROF)=XRECT+1.5-FLOAT(NRECT) RED(1)=RED(1)/2. SUM=0. DO 33 I=1,NPROF 33 SUM=SUM+RED(I) SUM=SUM*2. DO 34 I=1,NPROF RED(I)=RED(I)/SUM 34 BLUE(I)=RED(I) IF(X2.EQ.0.)GO TO 50 XRECT=FWHM2/2./VSTEP NRECT=XRECT+1.5 NPROF=NRECT DO 35 I=1,NPROF 35 RED2(I)=1. RED2(NPROF)=XRECT+1.5-FLOAT(NRECT) RED2(1)=RED2(1)/2. SUM=0. DO 36 I=1,NPROF 36 SUM=SUM+RED2(I) SUM=SUM*2. DO 37 I=1,NPROF RED2(I)=RED2(I)/SUM 37 BLUE2(I)=RED2(I) GO TO 50 C C INSTRUMENTAL PROFILE TABULATED AT SPECTRUM POINT SPACING. C RED AND BLUE WINGS BOTH START AT CENTRAL POINT. C THE PROFILE SHOULD SUM TO 1. 40 NPROF=X1 READ(5,41)(RED(I),I=1,NPROF) READ(5,41)(BLUE(I),I=1,NPROF) 41 FORMAT(10X,7F10.6) RED(1)=RED(1)/2. BLUE(1)=BLUE(1)/2. IF(X2.LE.0.)GO TO 50 READ(5,41)(RED2(I),I=1,NPROF) READ(5,41)(BLUE2(I),I=1,NPROF) RED2(1)=RED2(1)/2. BLUE2(1)=BLUE2(1)/2. 50 WRITE(6,51)(I,RED1(I),BLUE1(I),RED2(I),BLUE2(I),I=1,NPROF) 51 FORMAT(I5,4F10.6) WRITE(22)TEFF,GLOG,TITLE,WBEGIN,RESOLU,NWL,IFSURF,NMU,XMU,NEDGE, 1WLEDGE NH=(NWL+19999+19999)*NMU DO 52 I=1,NH 52 H(I)=0. WRITE(6,1117) 1117 FORMAT(1H1) IF(NMU.EQ.1)GO TO 150 DO 57 IWL=1,NWL READ(21)(QMU(I),I=1,NMU) DO 56 IMU=1,NMU Q=QMU(IMU) IWL1000=(IWL+20000)*NMU+IMU DO 53 I=1,NPROF 53 H(IWL1000-I*NMU)=H(IWL1000-I*NMU)+BLUE(I)*Q IWL998=(IWL+19998)*NMU+IMU DO 54 I=1,NPROF 54 H(IWL998+I*NMU)=H(IWL998+I*NMU)+RED(I)*Q 56 CONTINUE 57 CONTINUE IF(X2.EQ.0.)GO TO 160 REWIND 21 READ(21) DO 58 IWL=1,NWL IWLNMU=(IWL+19999)*NMU WT2=FLOAT(IWL-1)/FLOAT(NWL-1) WT1=1.-WT2 DO 58 IMU=1,NMU 58 H(IWLNMU+IMU)=H(IWLNMU+IMU)*WT1 DO 257 IWL=1,NWL READ(21)(QMU(I),I=1,NMU) WT2=FLOAT(IWL-1)/FLOAT(NWL-1) DO 256 IMU=1,NMU Q=QMU(IMU)*WT2 IWL1000=(IWL+20000)*NMU+IMU DO 253 I=1,NPROF 253 H(IWL1000-I*NMU)=H(IWL1000-I*NMU)+BLUE(I)*Q IWL998=(IWL+19998)*NMU+IMU DO 254 I=1,NPROF 254 H(IWL998+I*NMU)=H(IWL998+I*NMU)+RED(I)*Q 256 CONTINUE 257 CONTINUE GO TO 160 150 DO 157 IWL=1,NWL READ(21)QMU(1) Q=QMU(1) IWL1001=IWL+20001 DO 153 I=1,NPROF 153 H(IWL1001-I)=H(IWL1001-I)+BLUE(I)*Q IWL999=IWL+19999 DO 154 I=1,NPROF 154 H(IWL999+I)=H(IWL999+I)+RED(I)*Q 157 CONTINUE IF(X2.EQ.0.)GO TO 160 REWIND 21 READ(21) DO 358 IWL=1,NWL WT2=FLOAT(IWL-1)/FLOAT(NWL-1) WT1=1.-WT2 358 H(IWL+20000)=H(IWL+20000)*WT1 DO 357 IWL=1,NWL READ(21)QMU(1) WT2=FLOAT(IWL-1)/FLOAT(NWL-1) Q=QMU(1)*WT2 IWL1001=IWL+20001 DO 353 I=1,NPROF 353 H(IWL1001-I)=H(IWL1001-I)+BLUE(I)*Q IWL999=IWL+19999 DO 354 I=1,NPROF 354 H(IWL999+I)=H(IWL999+I)+RED(I)*Q 357 CONTINUE 160 REWIND 21 READ(21) DO 170 IWL=1,NWL READ(21)(QMU(IMU),IMU=1,NMU2) IWLNMU=(IWL+19999)*NMU DO 162 IMU=1,NMU 162 QMU(IMU)=H(IWLNMU+IMU) WRITE(22)(QMU(I),I=1,NMU2) IF(IWL.GT.LINOUT)GO TO 170 WAVE=WBEGIN*RATIO**(IWL-1) RESID=QMU(1)/QMU(NMU1) IRESID=RESID*1000.+.5 IPLOT=RESID*100.+1.5 IPLOT=MAX0(1,MIN0(101,IPLOT)) APLOT(IPLOT)='X' WRITE(6,2300)IWL,WAVE,IRESID,APLOT 2300 FORMAT(1H ,I5,F11.4,I7,101A1) APLOT(IPLOT)=' ' 170 CONTINUE READ(21)NLINES WRITE(22)NLINES DO 200 I=1,NLINES READ(21)LINDAT8,LINDAT4 WRITE(22)LINDAT8,LINDAT4 200 CONTINUE CALL EXIT END