SUBROUTINE INIBL0 C C AUXILIARY INITIALIZATION PROCEDURE C INCLUDE 'PARAMS.FOR' INCLUDE 'MODELP.FOR' INCLUDE 'LINDAT.FOR' INCLUDE 'SYNTHP.FOR' INCLUDE 'WINCOM.FOR' parameter (un=1.) character*2 iu character*6 ilab DIMENSION CROSS(MCROSS,MFRQ), * ABSO(MFREQ),EMIS(MFREQ),SCAT(MFREQ), * ABSOC(MFREQC),EMISC(MFREQC),SCATC(MFREQC) COMMON/LIMPAR/ALAM0,ALAM1,FRMIN,FRLAST,FRLI0,FRLIM COMMON/BLAPAR/RELOP,SPACE0,CUTOF0,TSTD,DSTD,ALAMC common/lasers/lasdel common/linrej/ilne(mdepth),ilvi(mdepth) common/velaux/velmax,iemoff,nltoff,itrad common/alsave/ALAM0s,ALASTs,CUTOF0s,CUTOFSs,RELOPs,SPACEs C C -------------------------------------------------------------- C Parameters controlling an evaluation of the synthetic spectrum C C -------------------------------------------------------------- C C ALAM0, ALAM1 - synthetic spectrum is evaluated between wavelengths C ALAM0 (initial) and ALAM1 (final), given in Anstroms C CUTOF0 - cutoff parameter for normal lines (given in Angstroms) C ie the maximum distance from the line center, in C which the opacity in the line is allowd to contribute C to the total opacity (recommended 5 - 10) C CUTOFS = SPACON C SPACON - spacing of the continuum wavelength points C (at the midpoint of teh total interval; actual spacing C is equidistant in log(lambda) C RELOP - the minimum value of the ratio (opacity in the line C center)/(opacity in continuum), for which is the line C taken into account (usually 1d-4 to 1d-3) C SPACE - the maximum distance of two neighbouring frequency C points for evaluating the spectrum; in Angstroms C C INLTE = 0 - pure LTE (no line in NLTE) C ne.0 - NLTE option, ie one or more lines treated C in the exact or approximate NLTE approach C IFHE2 gt.0 - He II line opacity in the first four series C (Lyman, Balmer, Paschen, Brackett) C for lines with lambda < 3900 A C is taken into account even if line list C does not contain any He II lines (i.e. C He II lines are treated as the hydrogen lines) C C IHYDPR = 0 - means that hydrogen lines Stark profiles C are calculated by approximate formulae C > 0 - hydrogen lines Stark profiles are calculated C in detail, using the Schoening & Butler tables; C (for 1-2 to 1-5; 2-3 to 2-10). C the tables are stored in file FOR0xx.dat, C where xx=IHYDPR; C higher Balmer lines are calculated as before C C the meaning of other parameters is quite analogous, for the C following lines C C IHE1PR - He I lines at 4471, 4026, 4387, and 4922 Angstroms C (tables calculated by Barnard, Cooper, and Shamey) C IHE2PR - for the He II lines calculated by Schoening and Butler, C if(ifeos.le.0) then READ(55,*) IFREQ,INLTE,ICONTL,INLIST,IFHE2 IF(LTE) INLTE=0 READ(55,*) IHYDPR,IHE1PR,IHE2PR READ(55,*) ALAM0,ALAST,CUTOF0,CUTOFS,RELOP,SPACE end if C IF(IDSTD.EQ.0) THEN ID1=5 NDSTEP=(ND-2*ID1)/2 IDSTD=2*ND/3 ELSE IF(IDSTD.LT.0) THEN ID1=1 NDSTEP=-IDSTD IDSTD=2*ND/3 END IF if(imode.le.-3) ndstep=1 c alam0s=alam0 alasts=alast cutof0s=cutof0 cutofss=cutofs relops=relop spaces=space C C if ALAST.lt.0 - set up vacuum wavelengths everywhere C vaclim=2000. if(alast.lt.0.) then alast=abs(alast) alasts=alast vaclim=1.e18 end if c if(inlte.lt.10) then lasdel=.true. else if(inlte.le.20) then inlte=inlte-10 lasdel=.false. else if(inlte.le.30) then inlte=inlte-20 ifreq=11 lasdel=.true. else if(inlte.le.40) then inlte=inlte-30 ifreq=11 lasdel=.false. end if C ibin(0)=mod(inlist,10) do ilist=1,mmlist tmlim(ilist)=tmolim ibin(ilist)=mod(inlist,10) ivdwli(ilist)=0 iun=19+ilist write(iu,622) iun 622 format(i2) amlist(ilist)='fort.' // iu end do c if(imode.ge.-3.and.imode.le.1) then nmlist=0 numlis=0 read(55,*,err=5,end=5) nmlist,(iunitm(ilist),ilist=1,nmlist) do ilist=1,nmlist write(iu,622) iunitm(ilist) amlist(ilist) ='fort.' // iu end do 5 continue c ilist=0 amlist(0)='fort.19' read(3,*,err=20,end=20) amlist(0),ibin(0) c ilist=0 10 continue ilist=ilist+1 read(3,*,end=20) amlist(ilist),ibin(ilist),tmlim(ilist) numlis=numlis+1 go to 10 20 continue if(numlis.gt.0) nmlist=numlis if(nmlist.gt.0.and.ifmol.eq.0) then write(*,*) 'NEEDS TO SET IFMOL > 0 with NMLIST>0' stop end if c ilist=0 ilab='ATOMIC' write(6,623) ilist,ilab,trim(amlist(ilist)),ibin(ilist) ilab='MOLEC ' do ilist=1,nmlist write(6,624) ilist,ilab,trim(amlist(ilist)),ibin(ilist), * tmlim(ilist) end do 623 format(/'************************'/ * ' LINE LISTS:'/ * /' ILIST',8x,'FILENAME IBIN TMLIM'/ * i4,2x,a6,2x,a,2x,i4,f11.1) 624 format( i4,2x,a6,2x,a,2x,i4,f11.1) end if c C c VTB - turbulent velocity (in km/s). In non-negative, this C value overwrites the value given by the standard input C read(55,*,err=30,end=30) VTB if(ifwin.le.0) then if(vtb.ge.0.) then WRITE(6,608) VTB 608 FORMAT(//' TURBULENT VELOCITY - CHANGED TO VTURB =', * 1PE10.3,' KM/S'/' ------------------'/) do id=1,nd vturb(id)=vtb*vtb*1.e10 end do end if end if C TSTD=TEMP(IDSTD) VTS=VTURB(IDSTD) DSTD=SQRT(1.4E7*TSTD+VTS) 30 continue C C angle points (in case the specific intensities are evaluated C C NMU0 - number of angles: C >0 - and if also ANG0>0, angles (mu's) equidistant C between 1 and ANG0 C >0 - and if also ANG0<0, angles (mu's) equidistant C between 0.7 and ANG0, and sinuses equidistatnt for C others C <0 - angles read in the next record C ANG0 - minimum mu (see above) C IFLUX - mode for evaluating angle-dependent intensities and C the corresponding flux: C =0 - no specifiec intensities are evaluated; only usual C flux is stored (unit 7 and 17) C =1 - specific intensities are evaluated; C and stored on unit 18 C =2 - (interesting only for the case of macroscopic C velocity field); specific intensities evaluated by C a simple formal solution (RESOLV) C NMU0=1 ANG0=1. ANGL(1)=1. WANGL(1)=0. IFLUX=0 velmax=3.e5 nltoff=0 iemoff=0 itrad=0 do id=1,nd wdil(id)=un end do if(ifwin.le.0) then READ(55,*,end=100,err=100) NMU0,ANG0,IFLUX C C determinantion of the angle points and weights C IF(NMU0.LT.0) THEN NMU0=IABS(NMU0) READ(55,*) (ANGL(IMU),IMU=1,NMU0) DO IMU=2,NMU0-1 WANGL(IMU)=0.5*(ANGL(IMU-1)+ANGL(IMU+1)) END DO WANGL(1)=0.5*(ANGL(1)-ANGL(2)) WANGL(NMU0)=0.5*(ANGL(NMU0-1)-ANGL(NMU0)) ELSE IF(ANG0.GT.0.) THEN IF(NMU0.GT.1) THEN DMU=(1.-ANG0)/(NMU0-1) DO IMU=1,NMU0 ANGL(IMU)=1.-(IMU-1)*DMU WANGL(IMU)=DMU END DO WANGL(1)=0.5*DMU WANGL(NMU0-1)=0.5*DMU WANGL(NMU0)=2.*DMU END IF ELSE ANGH=0.70710678 DMU=ANGH/(NMU0-1) DO IMU=1,NMU0 ANGL(IMU)=(IMU-1)*DMU ANGL(IMU)=SQRT(1.-ANGL(IMU)**2) IF(IMU.GT.1.AND.IMU.LT.NMU0) * WANGL(IMU)=0.5*(ANGL(IMU-1)+ANGL(IMU+1)) END DO WANGL(1)=0.5*(ANGL(1)-ANGL(2)) WANGL(NMU0)=0.5*(ANGL(NMU0-1)-ANGL(NMU0)) IF(ANG0.LT.0.) DMU=(ANGH+ANG0)/(NMU0-1) DO IMU=1,NMU0-2 ANGL(IMU+NMU0)=ANGH-IMU*DMU WANGL(IMU+NMU0)=DMU END DO WANGL(NMU0)=WANGL(NMU0)+0.5*DMU WANGL(2*NMU0-3)=0.5*DMU WANGL(2*NMU0-2)=2.*DMU NMU0=2*NMU0-2 END IF END IF IF(NMU0.LE.0) GO TO 100 WRITE(6,609) NMU0,(ANGL(I),I=1,NMU0) 609 FORMAT(//' SPECIFIC INTENSITIES COMPUTED FOR',I3, * ' ANGLES mu=cos(theta) ='/ * ' ---------------------------------', * '------------------------'// * (10F7.2)) 100 CONTINUE else itrad=1 read(55,*,end=110,err=110) velmax,ITRAD,nltoff,iemoff 110 write(6,602) velmax,itrad,nltoff,iemoff if(velmax.lt.0.) then velmax=3.e5 go to 120 end if 602 format(//' velmax (velocity for line rejection)', * ' itrad,nltoff,iemoff',f10.1,2i3) C C Set up rays and weights C call velset call radtem CALL SETRAY CALL WGTJH1 C end if C 120 CONTINUE velmax=velmax*1.e5 do id=1,nd ilvi(id)=0 ilne(id)=0 if(vel(id).gt.velmax.and.iemoff.eq.0) ilvi(id)=1 if(vel(id).gt.velmax.and.nltoff.gt.0.and.iemoff.gt.0) * ilne(id)=1 end do C IF(IMODE.EQ.-1) THEN INLTE=0 CUTOF0=0. END IF C C continuum frequencies C if(ifwin.le.0) then alam0=alam0s if(alam0s.eq.0.) alam0=5.e7/temp(1)/10. if(alam0s.lt.0.) alam0=-5.e7/temp(1)/alam0s alast=alasts if(alasts.eq.0.) alast=5.e7/temp(1)*20. if(alasts.lt.0.) alast=-5.e7/temp(1)*alasts c if(alast.gt.1.e5) alast=1.e5 ALAMC=(ALAM0+ALAST)*0.5 if(space.eq.0.) space=4.3e-8*sqrt(temp(idstd))*alamc if(space.lt.0.) space=-5.72e-8*sqrt(temp(idstd))*alamc*space SPACF=2.997925E18/ALAMC/ALAMC*SPACE WRITE(6,601) ALAM0,ALAST,CUTOF0,RELOP,SPACF,SPACE CUTOF0=0.1*CUTOF0 SPACE0=SPACE*0.1 ALAM0=1.D-1*ALAM0 ALAST=1.D-1*ALAST ALAMC=ALAMC*0.1 ALST00=ALAST FRLAST=2.997925D17/ALAST NFREQ=2 FREQ(1)=2.997925D17/ALAM0 FREQ(2)=FRLAST C else C spacon=cutofs IF(SPACON.EQ.0) SPACON=3. XFR=(ALAST-ALAM0)/SPACON NFREQC=int(XFR)+1 NFREQC=MIN(NFREQC,MFREQC) NFREQC=MAX(NFREQC,2) DLAMLO=LOG10(ALAST/ALAM0)/(NFREQC-1) AL0L=LOG10(ALAM0) alambe=alam0 DO IJ=1,NFREQC AL=AL0L+(IJ-1)*DLAMLO ALAM=EXP(2.3025851*AL) WLAMC(IJ)=ALAM FREQC(IJ)=2.997925E18/ALAM END DO ALAMC=(ALAM0+ALAST)*0.5 SPACF=2.997925E18/ALAMC/ALAMC*SPACE WRITE(6,601) ALAM0,ALAST,CUTOF0,RELOP,SPACF,SPACE CUTOF0=0.1*CUTOF0 SPACE0=SPACE*0.1 ALAM0=1.D-1*ALAM0 ALAST=1.D-1*ALAST ALAMC=ALAMC*0.1 ALST00=ALAST FRLAST=2.997925D17/ALAST NFREQ=2 FREQ(1)=2.997925D17/ALAM0 FREQ(2)=FRLAST c end if c CALL SIGAVS IF(IHYDPR.NE.0) THEN CALL HYDINI CALL XENINI END IF IF(IHE1PR.GT.0) CALL HE1INI IF(IHE2PR.GT.0) CALL HE2INI C C auxiliary quantities for dissolved fractions C DO ID=1,ND CALL DWNFR0(ID) CALL WNSTOR(ID) END DO C c pretabulate expansion coefficients for the Voigt function c CALL PRETAB c c calculate the characteristic standard opacity c IF(IMODE.LE.2) THEN if(ifwin.le.0.and.ndstep.eq.0) then c c old procedure c CALL CROSET(CROSS) DO ID=1,ND CALL OPAC(ID,CROSS,ABSO,EMIS,SCAT) ABSTD(ID)=MIN(ABSO(1),ABSO(2)) END DO else c c new procedure c if(ifwin.le.0) then nfreqc=ifix(real(cutofs,4)) if(nfreqc.eq.0) nfreqc=mfreq all0=log(alam0) all1=log(alast) dlc=(all1-all0)/(nfreqc-1) do ijc=1,nfreqc wlamc(ijc)=exp(all0+(ijc-1)*dlc) freqc(ijc)=2.997925e17/wlamc(ijc) end do CALL CROSEW(CROSS) do id=1,nd CALL OPACON(ID,CROSS,ABSOC,EMISC,SCATC) do ijc=1,nfreqc abstdw(ijc,id)=absoc(ijc) end do end do c write(*,*) 'abstdw(1,ij)',(abstdw(ij,1),ij=1,nfreqc) c write(*,*) 'abstdw(50,ij)',(abstdw(ij,50),ij=1,nfreqc) c else CALL CROSEW(CROSS) DO ID=1,ND CALL OPACW(ID,CROSS,ABSO,EMIS,ABSOC,EMISC,SCATC,0) DO IJ=1,NFREQC ABSTDW(IJ,ID)=ABSOC(IJ)/DENSCON(ID) END DO END DO end if end if END IF C 601 FORMAT(//'----------------------------------------------'/ * ' BASIC INPUT PARAMETERS FOR SYNTHETIC SPECTRA'/ * ' ---------------------------------------------'/ * ' INITIAL LAMBDA',28X,1H=,F10.3,' ANGSTROMS'/ * ' FINAL LAMBDA',28X,1H=,F10.3,' ANGSTROMS'/ * ' CUTOFF PARAMETER',26X,1H=,F10.3,' ANGSTROMS'/ * ' MINIMUM VALUE OF (LINE OPAC.)/(CONT.OPAC) =',1PE10.1/ * ' MAXIMUM FREQUENCY SPACING',17X,1H=,1PE10.3,' I.E.', * 0PF6.3,' ANGSTROMS'/ * ' ---------------------------------------------'/) c write(6,612) idstd,ndstep 612 format(/'IDSTD, NDSTEP = ',2i5/) RETURN END