Module arguments ! Implicit None ! Character(len=256),dimension(5) :: c_arg = "" ! Character(len=5) :: EltIon = "" Real :: lambda = Huge(0E0) ! Angstroem Real :: Teff = Huge(0E0) ! Kelvin Real :: logg = Huge(0E0) ! log(cm/s2) Real :: Fe_H = Huge(0E0) ! logarithmic scale with respect to the solar value and where H is 12. ! Integer :: lambda_ind_i, lambda_ind_j ! Integer, parameter :: Nelt=3 ! Number of elements in given of ionization stage Character(len=5),dimension(Nelt),parameter :: Data_EltIon=(/"MgI ","CaI ","CaII "/) Real,dimension(Nelt,3),parameter :: Data_lambda=Reshape((/8473.688,8715.255,8736.012,& 8525.707,8583.318,8633.929,& 8498.018,8542.086,8662.135/),(/Nelt,3/),order=(/2,1/)) Real,parameter :: Tol_lambda = 1.0 Real,parameter :: Teffmin = 3500, Teffmax = 5500 Real,parameter :: loggmin = 0.5, loggmax = 2.0 Real,parameter :: Fe_Hmax = 0.5 Real,dimension(Nelt,3),parameter :: Fe_Hmin = Reshape((/-1.0,-1.5,-2.0, -1.0,-1.0,-1.0, -4.0,-4.0,-4.0/),(/3,3/),order=(/2,1/)) End Module arguments ! !*********************************************************************** ! Module coeff_fits ! USE arguments, ONLY: Nelt ! Implicit None ! Real,dimension(Nelt,3,20), parameter :: mat_coeff=Reshape((/& 1.03E+00,3.32E-02,-1.79E-02,-1.48E-01,1.50E-02,0.,0.,2.24E-02,4.39E-02,0.,-2.31E-02,-6.92E-03,0.,0.,4.43E-02,-1.23E-02,0., & ! MgI 8473 -2.26E-02,1.90E-02,3.75E-02, & ! MgI 8473 (END) 1.22E+00,1.43E-01,-9.26E-02,-3.90E-02,-1.04E-01,6.79E-02,-2.76E-02,-1.86E-02,9.46E-02,-2.95E-01,0.,-2.49E-02,0.,0., & ! MgI 8715 1.84E-01,-6.67E-02,2.82E-02,-7.15E-02,0.,1.72E-01, & ! MgI 8715 (END) 1.31E+00,-6.84E-02,-2.36E-02,1.61E-01,-6.34E-02,0.,0.,8.05E-02,-2.16E-02,-2.42E-01,3.44E-02,-2.50E-02,0.,0.,1.20E-01,0., & ! MgI 8736 0.,-1.18E-01,5.08E-02,-5.34E-02, & ! MgI 8736 (END) 5.95E-01,-2.81E-02,-1.67E-02,7.11E-02,4.74E-02,2.41E-02,1.22E-02,2.09E-01,5.01E-02,-2.49E-02,-9.98E-02,1.61E-02,0.,0.,0., & ! CaI 8525 0.,0.,1.04E-02,1.06E-02,6.04E-02, & ! CaI 8525 (END) 7.05E-01,0.,0.,2.76E-02,5.13E-02,2.12E-02,0.,1.84E-01,3.87E-02,-3.69E-02,-1.02E-01,0.,0.,0.,0.,0.,0.,1.14E-02,6.51E-03, & ! CaI 8583 6.25E-02, & ! CaI 8583 (END) 7.97E-01,6.00E-02,0.,0.,6.63E-02,2.02E-02,0.,1.32E-01,2.81E-02,-3.44E-02,-1.20E-01,0.,0.,0.,0.,0.,0.,1.53E-02,0.,4.85E-02, & ! CaI 8633 0.988,0.0368,0.00915,-0.0844,0.0422,0.,0.,-0.0387,-0.0381,0.112,-0.0254,0.0115,0.,0.,-0.0322,0.,0.,0.0417,0.0293,-0.0747, & ! CaII 8498 A 0.972,0.0371,0.00305,-0.0274,0.0466,0.,0.00645,0.00533,-0.0293,0.0777,-0.0243,0.0142,-0.0108,0.,-0.0532,0.,-0.00581,0., & ! CaII 8542 A 0.0295,-0.0774, & ! CaII 8542 A (END) 0.973,0.0382,0.00221,-0.0318,0.0460,0.,0.,-0.0167,-0.0291,0.0938,-0.0232,0.0134,-0.00903,0.,-0.0547,0.,0.,0.0207,0.0321, & ! CaII 8662 A -0.0814/), & ! CaII 8662 A (END) (/Nelt,3,20/),order=(/3,2,1/)) End Module coeff_fits ! !*********************************************************************** ! Program NLTE_Gaia_lines ! ! Give the NLTE/LTE Equivalent Width (EW), W/W*(line), interpolated at ! Teff, logg and [Fe/H] for MgI, CaI and CaII lines in the RVS/Gaia, ! using the coefficient fits given in Merle et al. (2011) MNRAS ! The program needs 5 arguments: ! element with ionization stage, lambda (A), Teff (K), log g (log cgs) and [Fe/H] ! If: $ gfortran -o NLTE_Gaia_lines NLTE_Gaia_lines.f90 ! Then: $ ./NLTE_Gaia_lines CaII 8498 5200 1.35 -0.5 ! ! This means that we want to compute the NLTE/LTE EW for the CaII 8498.018 A line ! for the model atmosphere with Teff = 5200 K, log g = 1.35 and [Fe/H] = -0.5 ! ! 0.0 when [Fe/H] > 0.0 ! We remind that [alpha/Fe] = 0.4[Fe/H] when -1.0 < [Fe/H] < 0.0 ! 0.4 when [Fe/H] < -1.0 ! !History: 2011-07-28 (TM) Creation ! USE arguments USE coeff_fits ! Implicit None ! Integer :: Narg,ios Real :: WWS=Huge(0E0) ! Narg = Command_argument_count() ! Call CMD(Narg) Call VAL(EltIon,lambda,Teff,logg,Fe_H) Call FIT(EltIon,lambda,Teff,logg,Fe_H,WWS) ! !Write(*,"(2A,F6.0,A,F5.2)")"WWS(",Trim(EltIon),lambda,") = ", WWS Write(*,"(F5.2)") WWS ! CONTAINS ! !-------------------------------------------------------------------- ! Real Function POLY_XYZ(x,y,z,coeff) Real,intent(in) :: x,y,z Real,dimension(:),intent(in) :: coeff ! POLY_XYZ = coeff(1) + coeff(2)*x + coeff(3)*y + coeff(4)*z + & coeff(5)*x**2 + coeff(6)*x*y + coeff(7)*y**2 + coeff(8)*x*z + & coeff(9)*y*z + coeff(10)*z**2 + coeff(11)*x**3 + coeff(12)*x**2*y + & coeff(13)*x*y**2 + coeff(14)*y**3 + coeff(15)*x**2*z + coeff(16)*x*y*z + & coeff(17)*y**2*z + coeff(18)*x*z**2 + coeff(19)*y*z**2 + coeff(20)*z**3 End Function POLY_XYZ ! !-------------------------------------------------------------------- ! Subroutine FIT(EltIon,lambda,Teff,logg,Fe_H,WWS) Character(len=*),intent(in) :: EltIon Real,intent(in) :: lambda,Teff,logg,Fe_H Real,intent(out) :: WWS ! Integer :: I Real :: x,y,z, xmid,ymid,zmid, xdev,ydev,zdev Real,dimension(20) :: coeff Logical :: blabla=.FALSE. ! xmid = 0.5*(Teffmax+Teffmin) xdev = 0.5*(Teffmax-Teffmin) ymid = 0.5*(loggmax+loggmin) ydev = 0.5*(loggmax-loggmin) zmid = 0.5*(Fe_Hmax+Fe_Hmin(lambda_ind_i,lambda_ind_j)) zdev = 0.5*(Fe_Hmax-Fe_Hmin(lambda_ind_i,lambda_ind_j)) ! IF (blabla) Write(*,*)"xmid,xdev,ymid,ydev,zmid,zdev:",xmid,xdev,ymid,ydev,zmid,zdev ! x = (Teff-xmid)/xdev y = (logg-ymid)/ydev z = (Fe_H-zmid)/zdev ! IF (blabla) Write(*,*)"x,y,z:",x,y,z ! coeff = mat_coeff(lambda_ind_i,lambda_ind_j,:) ! IF (blabla) Write(*,'(ES10.2)') (coeff(I),I=1,20) ! WWS = POLY_XYZ(x,y,z,coeff) ! End Subroutine FIT ! !-------------------------------------------------------------------- ! Subroutine CMD(Narg) Integer,intent(in) :: Narg Integer :: I Logical :: blabla = .FALSE. Character(Len=10) :: EltIon_line ! IF (Narg/=4) Then Call INFO() STOP "NLTE_GAIA_LINES PROGRAM: PB-00: MISSING ARGUMENTS." End If ! Do I=1,Narg Call Get_command_argument(I,c_arg(I),status=ios) IF (ios/=0) STOP "NLTE_GAIA_LINES PROGRAM: PB-01: FAILED ARGUMENTS." End Do ! IF(blabla) Write(*, *) "c_arg:", ( (Trim(c_arg(I))//" "),I=1,5 ) ! Read(c_arg(1),"(A)",iostat=ios) EltIon_line IF (ios/=0) STOP "NLTE_GAIA_LINES PROGRAM: PB-02: ELT WITH IONIZATION STAGE, ONLY MgI, CaI and CaII AVAILABLE." Select Case (EltIon_line) Case ("MgI_8473") EltIon = "MgI " ; lambda = 8473.0 Case ( "MgI_8715 " , "MgI_8710 " , "MgI_8712 " , "MgI_8717 " ) EltIon = "MgI " ; lambda = 8715.0 Case ("MgI_8736") EltIon = "MgI " ; lambda = 8736.0 Case ("CaI_8525") EltIon = "CaI " ; lambda = 8525.0 Case ("CaI_8583") EltIon = "CaI " ; lambda = 8583.0 Case ("CaI_8633") EltIon = "CaI " ; lambda = 8633.0 Case ("CaII_8498 ") EltIon = "CaII " ; lambda = 8498.0 Case ("CaII_8542 ") EltIon = "CaII " ; lambda = 8542.0 Case ("CaII_8662") EltIon = "CaII " ; lambda = 8662.0 Case Default Write(*,*) " NLTE_GAIA_LINES PROGRAM: PB-03: WRONG ARGUMENT" Write(*,*) " EltIon_line AVAILABLE ARE : " Write(*,*) " MgI_8473 , MgI_8715 (MgI_8710, MgI_8712, MgI_8717), MgI_8736 ," Write(*,*) " CaI_8525 , CaI_8583 , CaI_8633 ," Write(*,*) " CaII_8498 , CaII_8542 , CaII_8662." STOP End Select ! Read(c_arg(2),"(F10.0)",iostat=ios) Teff IF (ios/=0) STOP "NLTE_GAIA_LINES PROGRAM: PB-04: EFFECTIVE TEMPERATURE IN KELVIN IS REQUIRED." Read(c_arg(3),"(F5.0)",iostat=ios) logg IF (ios/=0) STOP "NLTE_GAIA_LINES PROGRAM: PB-05: LOG OF SURFACE GRAVITY IN LOG(CGS) REQUIRED." Read(c_arg(4),"(F5.0)",iostat=ios) Fe_H IF (ios/=0) STOP "NLTE_GAIA_LINES PROGRAM: PB-11: GLOBAL METALLICITY WITH RESPECT TO THE SUN IS REQUIRED." ! IF(blabla) Write(*,*) EltIon,lambda,Teff,logg,Fe_H End Subroutine CMD ! !-------------------------------------------------------------------- ! Subroutine VAL(EltIon,lambda,Teff,logg,Fe_H) Character(len=*),intent(in) :: EltIon Real,intent(in) :: lambda,Teff,logg,Fe_H ! Logical,dimension(5) :: larg=.FALSE. Logical,parameter :: blabla=.FALSE. Integer :: I,J ! larg(1)=Any(Data_EltIon==Trim(EltIon)) ! If (larg(1)) Then IF(blabla) Write(*,*)"Elt: ",EltIon," OK." ! Do I=1,Size(Data_EltIon) IF (EltIon==Data_EltIon(I)) EXIT End Do ! lambda_ind_i = I ! larg(2) = Any(Abs(Data_lambda(I,:)-lambda)<=Tol_lambda) J=Minloc(Abs(Data_lambda(I,:)-lambda),dim=1) lambda_ind_j = J ! If (EltIon=="MgI" .AND. (Abs(lambda-8710.169)<= Tol_lambda .OR. & Abs(lambda-8712.678)<= Tol_lambda .OR. & Abs(lambda-8717.812)<= Tol_lambda)) Then larg(2) = .TRUE. lambda_ind_j = 2 End If ! If (.NOT.larg(2)) STOP " " IF (blabla) Write(*,*)"Line: ",lambda," A, OK." larg(3) = Teff>=Teffmin .AND. Teff<=Teffmax If (.NOT.larg(3)) STOP "NLTE_GAIA_LINES PROGRAM: PB-08: TEFF OUTSIDE OF THE VALID RANGE: 3500 < TEFF < 5500." IF (blabla) Write(*,*)"Teff = ",Teff," K, OK." larg(4) = logg>=loggmin .AND. logg<=loggmax If (.NOT.larg(4)) STOP "NLTE_GAIA_LINES PROGRAM: PB-09: LOGG OUTSIDE OF THE VALID RANGE: 0.5 < LOGG < 2.0." IF (blabla) Write(*,*)"log g = ",logg," OK." ! larg(5) = Fe_H>=Fe_Hmin(I,J) .AND. Fe_H<=Fe_Hmax If (.NOT.larg(5)) STOP "NLTE_GAIA_LINES PROGRAM: PB-10: [Fe/H] OUTSIDE OF THE VALID RANGE: [Fe/H]min < [Fe/H] < 0.5." IF (blabla) Write(*,*)"[Fe/H] = ",Fe_H," OK." ! Else STOP "NLTE_GAIA_LINES PROGRAM: PB-06: NO FITS FOR THIS ELEMENT." End If ! End Subroutine VAL ! !-------------------------------------------------------------------- ! Subroutine INFO() Write(*,*)"" Write(*,*)" This fortran program give the NLTE/LTE Equivalent Width (EW), W/W*(line)," Write(*,*)" interpolated at Teff, logg and [Fe/H] for MgI, CaI and CaII lines in the RVS/Gaia," Write(*,*)" using the coefficient fits given in Merle et al. (2011) MNRAS" Write(*,*)" The program needs 4 arguments:" Write(*,*)" element with ionization stage and lambda (A), Teff (K), log g (log cgs) and [Fe/H]." Write(*,*) Write(*,*)" If: $ gfortran -o NLTE_Gaia_lines NLTE_Gaia_lines.f90" Write(*,*)" Then: $ ./NLTE_Gaia_lines CaII_8498 5200 1.35 -0.5" Write(*,*) Write(*,*)" This means that you want to compute the NLTE/LTE EW for the CaII 8498.018 A line" Write(*,*)" for the model atmosphere with Teff = 5200 K, log g = 1.35 and [Fe/H] = -0.5" Write(*,*) Write(*,*)" Elements and lines available:" Write(*,*)" MgI: 8473, <8715> (8710, 8712, 8717), 8736" Write(*,*)" CaI: 8525, 8583 , 8633" Write(*,*)" CaII: 8498, 8542 , 8662" Write(*,*) Write(*,*)" Thibault Merle -- 2011 07 28 --" Write(*,*) End Subroutine INFO ! End Program NLTE_Gaia_lines