Commit 9d338631 authored by Johann Cohen-Tanugi's avatar Johann Cohen-Tanugi
Browse files

Merge branch '41-improvements-to-sed-class' into u/jct/python

parents 324e4a30 d2b8d46c
import os
from ._lephare import *
from ._flt import *
from .filter import *
from .sedtolib import *
from .mag_gal import *
from .zphota import *
from .filter import *
from .sedtolib import *
......
......@@ -54,7 +54,7 @@ class MagGal(Runner):
return
if __name__ == "__main__":
if __name__ == "__main__":
runner = MagGal()
runner.run()
runner.end()
This diff is collapsed.
This diff is collapsed.
......@@ -34,15 +34,15 @@ public:
vector<double> mag,flux, ///< magnitude and flux of the model
flux_ori, ///< original flux vector before any transformation is applied
mag_ori; ///< original mag vector before any transformation is applied
string name,em;
int jtype,nummod;
int rec,nlib,index,index_z0,nbusIR;
static int nrec;
string name;
bool has_emlines; ///< True if the emission lines have been computed, false if not
int nummod;
int nlib,index,index_z0,nbusIR;
double red,redn,redp,chi2,prob,dm,lnir,luv,lopt;
double mass,age,sfr,ltir; // need to put it out of GalSED since used in the PDF without knowing that it's a gal.
double ebv,mag0;
int extlawId;
double qi[4]; // Cedric #photons.cm-2/s ionisant pour 4 elements
double qi[4]; ///< Store the number flux (phot/cm\f$^{-2}\f$s\f$^{-1}\f$) of ionizing photons for HeII, HeI, H, and H2. See SED::calc_ph. In practice, qi[2] only is used, and only for the physical modeling of emission lines (EM_LINES="PHYS", see GalMag::read_SED)
vector <oneElLambda> fac_line,fac_line_ori;
// Constructors defined in SED.cpp
......@@ -59,14 +59,11 @@ public:
flux_ori = p.flux_ori;
mag_ori = p.mag_ori;
name = p.name;
em = p.em;
jtype = p.jtype;
has_emlines = p.has_emlines;
nummod = p.nummod;
rec = p.rec;
nlib = p.nlib;
index = p.index;
index_z0 = p.index_z0;
nrec = p.nrec;
red = p.red;
redn = p.redn;
redp = p.redp;
......@@ -128,16 +125,51 @@ public:
virtual void fit(const onesource & source,const int imagm, const string zfix, const double consideredZ, const double lmasi, const double lmass,const double magabsB[2], const double magabsF[2], const double funz0, const int bp[2]){};
virtual void fitIR(const onesource& source,const int imagm, const string zfix, const double consideredZ, const string fit_frsc, cosmo lcdm){};
virtual void writeSED(ofstream& ofs, ofstream& ofsPhys, ofstream& ofsDoc);
virtual void writeLib(string outasc, ofstream& ofsBin, ofstream& ofsDat, vector<flt> allFilters, string magtyp){};
inline void writeSED(const string& binFile, const string& physFile, const string& docFile) {
ofstream sdocOut, sphysOut, sbinOut;
sdocOut.open(docFile.c_str());
if(!sdocOut){
throw invalid_argument("Can't open doc file " + docFile);
}
sbinOut.open(binFile.c_str(), ios::binary | ios::out);
if(!sbinOut){
throw invalid_argument("Can't open bin file " + binFile) ;
}
if(type[0]=='G' || type[0]=='g'){
sphysOut.open(physFile.c_str());
if(!sphysOut){
throw invalid_argument("Can't open phys file " + physFile);
}
}
writeSED(sbinOut, sphysOut, sdocOut);
};
virtual void writeMag(bool outasc, ofstream& ofsBin, ofstream& ofsDat, vector<flt> allFilters, string magtyp){};
inline void readSEDBin(const string& fname) {
ifstream sbinIn;
sbinIn.open(fname.c_str(), ios::binary);
if(!sbinIn){throw invalid_argument("Can't open file " + fname); }
readSEDBin(sbinIn);
}
/// read the SED library when it is in binary format
virtual void readSEDBin(ifstream& ins);
virtual void readLib(ifstream& ins){};
virtual void readMagBin(ifstream& ins){};
virtual void sumEmLines(){};
/// for each magnitude \a #mag[k] compute kcorr = mag[k] - mag_z0[k] - distMod
virtual void kcorrec(){};
virtual void prop(){};
virtual void add_neb_cont(){}; // Add continuum
/*!
* Compute the number flux of photons able to ionize HeII, HeI, H, and H2
* For a given SED, this amounts to compute the integral
* \f$\int_0^{w_i} SED(\lambda)\cdot \frac{\lambda}{hc}\,d\lambda\quad,\f$
* where \f$w_i\f$=54.42, 24.52, 13.60, and 1108.7 A for HeII, HeI, H, and H2 respectively,
* and where \f$hc\f$ is in ergs.A. This normalization assumes that the SED are provided in args/cm2/s/A.
* In practice the integral is approximated by :
\f$\sum_{\lambda_{min}}^{w_k}\frac{SED_{j-1}+SED_j}{2}\cdot(\lambda_j-\lambda_{j-1})\cdot\frac{\lambda_j}{hc}\f$.
*
* Results are stored in the q_i array member of size 4 of the SED instance.
*/
virtual void calc_ph(){}; // Number of inoizing photons
virtual void SEDproperties(){};
......@@ -149,6 +181,7 @@ public:
void fit_chi2(const onesource& source, const int imagm);
int check_discrepent(const onesource& source, const int imagm,const double levelChi2,const double levelMag);
void fit_chi2IR(const onesource& source, const int imagm, const string fit_frsc, const double dmcor);
inline bool is_same_model(const SED& other) {return ((*this).nummod == other.nummod && (*this).ebv == other.ebv) ;}
};
......@@ -200,7 +233,7 @@ public:
};
// Constructors defined in SED.cpp
GalSED(const string nameC, int nummodC=0, string typeC="N");
GalSED(const string nameC, int nummodC=0);
GalSED(const string nameC, double tauC, double ageC, int keepAgeC, string formatC, int nummodC, string typeC, int idAgeC);
~GalSED();
......@@ -211,7 +244,6 @@ public:
void generateEmPhys(double zmet, double qi);
void generateEmSpectra();
void generateEmSpectra(int nstep);
void fillPhys(const vector<double>& phys);
void sumEmLines();
void kcorrec();
void rescaleEmLines();
......@@ -222,8 +254,8 @@ public:
void writeSED(ofstream& ofs, ofstream& ofsPhys, ofstream& ofsDoc);
void readSEDBin(ifstream& ins);
void writeLib(string outasc, ofstream& ofsBin, ofstream& ofsDat, vector<flt> allFilters, string magtyp);
void readLib(ifstream& ins);
void writeMag(bool outasc, ofstream& ofsBin, ofstream& ofsDat, vector<flt> allFilters, string magtyp);
void readMagBin(ifstream& ins);
void generate_spectra(double zin, double dmin, vector<opa> opaAll);
void fit(const onesource& source,const int imagm, const string zfix, const double consideredZ, const double lmasi, const double lmass,const double magabsB[2], const double magabsF[2], const double funz0, const int bp[2]);
......@@ -240,13 +272,13 @@ class QSOSED : public ExtSED{
public:
QSOSED(SED const& p):ExtSED(p) {};
QSOSED(QSOSED const& p):ExtSED(p) {};
QSOSED(const string nameC, int nummodC=0, string typeC="N");
QSOSED(const string nameC, int nummodC=0);
~QSOSED();
void writeLib(string outasc, ofstream& ofsBin, ofstream& ofsDat, vector<flt> allFilters, string magtyp);
void writeMag(bool outasc, ofstream& ofsBin, ofstream& ofsDat, vector<flt> allFilters, string magtyp);
void clean(){ExtSED::clean();};
void readLib(ifstream& ins);
void readMagBin(ifstream& ins);
void prop();
void generate_spectra(double zin, double dmin, vector<opa> opaAll);
......@@ -262,11 +294,11 @@ class StarSED : public SED{
public:
StarSED(SED const& p):SED(p) {};
StarSED(StarSED const& p):SED(p) {};
StarSED(const string nameC, int nummodC=0, string typeC="N");
StarSED(const string nameC, int nummodC=0);
~StarSED();
void writeLib(string outasc, ofstream& ofsBin, ofstream& ofsDat, vector<flt> allFilters, string magtyp);
void readLib(ifstream& ins);
void writeMag(bool outasc, ofstream& ofsBin, ofstream& ofsDat, vector<flt> allFilters, string magtyp);
void readMagBin(ifstream& ins);
void prop();
void generate_spectra(double dmin);
......
......@@ -50,7 +50,6 @@ public:
///write time of creation and the number of SED recorded to the doc file
void print_time_tofile(time_t result){
sdocOut << "CREATION_DATE " << asctime(std::localtime(&result));
sdocOut << "NUMBER_ROWS " << (SED::nrec) << endl;
}
/// print config and type onscreen and in the doc file
virtual void print_info();
......@@ -131,7 +130,7 @@ void SEDLib<T>::write_SED_lib(){
template <class T>
void SEDLib<T>::readSED(string sedFile, string sedFormat, int nummod, string type){
// Create one object "SED" and fill it with one ascii file
T oneSEDascii(sedFile,nummod,type);
T oneSEDascii(sedFile,nummod);
oneSEDascii.read(sedFile);
//The vector has only one element SED since only one age
allSED.push_back(oneSEDascii);
......@@ -306,7 +305,7 @@ void SEDLib<GalSED>::readSED(string sedFile, string sedFormat, int nummod, strin
}else{
//ASCII by default
GalSED oneSEDascii(sedFile,nummod,type);
GalSED oneSEDascii(sedFile,nummod);
oneSEDascii.read(sedFile);
resultSED.push_back(oneSEDascii);
}
......
......@@ -129,7 +129,6 @@ PYBIND11_MODULE(_lephare, mod) {
/******** CLASS SED *********/
py::class_<SED>(mod, "SED")
.def_readonly_static("nrec", &SED::nrec)
.def(py::init<const string, int, string>(), py::arg("name"), py::arg("nummod") = 0, py::arg("type") = "N")
.def(py::init<const string, double, double, int, string, int, string, int>(), py::arg("name"), py::arg("tau"), py::arg("age"), py::arg("keepAge"), py::arg("format"), py::arg("nummod"), py::arg("type"), py::arg("idAge"))
.def(py::init<const SED>())
......@@ -145,11 +144,26 @@ PYBIND11_MODULE(_lephare, mod) {
.def("writeType", &SED::writeType)
.def("zaround", &SED::zaround)
.def("indexz", &SED::indexz)
.def("sumSpectra", &SED::read)
.def("readSEDBin", static_cast<void (SED::*)(const string&)>(&SED::readSEDBin))
.def("writeSED", static_cast<void (SED::*)(const string&, const string&, const string&)>(&SED::writeSED))
.def("data",
[](const SED &f){
int N = f.lamb_flux.size();
py::array_t<double> result = py::array_t<double>(N*2);
py::buffer_info buf = result.request();
double *ptr = (double *) buf.ptr;
for(size_t i = 0; i < N; i++){
ptr[i] = f.lamb_flux[i].lamb;
ptr[N+i] = f.lamb_flux[i].val;
}
result.resize({N,N});
return result;
}
)
;
py::class_<StarSED, SED>(mod, "StarSED")
.def(py::init<const string, int, string>(), py::arg("name"), py::arg("nummodC") = 0, py::arg("typeC") = "S")
.def(py::init<const string, int>(), py::arg("name"), py::arg("nummodC") = 0)
.def("prop", &StarSED::prop)
.def("read", &StarSED::read)
;
......
......@@ -65,7 +65,7 @@ Mag::Mag(keymap & key_analysed){
throw runtime_error("You are probably using the old parametrisation of Z_STEP since Z MIN > Z MAX in Z_STEP. Stop here. ");
}
// Output file in ascii ?
outasc = ((key_analysed["LIB_ASCII"]).split_string("NO",1))[0];
outasc = ((key_analysed["LIB_ASCII"]).split_bool("NO",1))[0];
// Want to display the template number on the screen
// VERBOSE output file - YES default
......@@ -130,7 +130,7 @@ void Mag::open_files(){
if(!sopaListIn){throw invalid_argument("Can't open file " + opaListFile); }
// Open an ascii file to store the predicted magnitudes
if(outasc[0] == 'Y' ||outasc[0] == 'y'){
if(outasc){
// Name of the output ascii file
datFile=lepharework + "/lib_mag/" + colib + ".dat";
sdatOut.open(datFile.c_str());
......@@ -281,7 +281,6 @@ void Mag::write_doc(){
sdocOut << "LIB_TYPE STAR" << endl;
}
sdocOut << "LIB_NAME " << lib << endl;
sdocOut << "NUMBER_ROWS " << (SED::nrec) << endl;
sdocOut << "FILTER_FILE " << filtFile << endl;
sdocOut << "FILTERS " ;
for (vector<flt>::iterator itf=allFlt.begin(); itf<allFlt.end(); ++itf)sdocOut << itf->name << "," ;
......@@ -300,7 +299,7 @@ void Mag::write_doc(){
sdocOut << endl << "EB_V " ;
for (int k=0 ; k<nebv ; k++){sdocOut << ebv[k] << "," ;};
sdocOut << endl << "EM_LINES " << emlines << endl;
sdocOut << "LIB_ASCII " << outasc << endl;
sdocOut << "LIB_ASCII " << (outasc?"YES":"NO") << endl;
time_t result = time(nullptr);
sdocOut << "CREATION_DATE " << asctime(std::localtime(&result));
......@@ -365,8 +364,6 @@ void GalMag::close_files(){
// Create the emission lines
// Compute the modeled magnitude. Done here to optimize the library size
void GalMag::read_SED(){
int nrec,jtype=0;
double age,luv,lopt,lnir,ltir,mass, sfr, zmet, tau, d4000,qi,MNUV_int,nuvr,val;
// double luv_int,luv_ext,kEBV;
GalSED saveSED("bid.dat");
......@@ -402,7 +399,7 @@ void GalMag::read_SED(){
ssedDocIn.ignore(1000000,'\n');
// read the physical parameters linked to this SED.
ssedPhysIn >> nrec >> jtype >> age >> luv >> lopt >> lnir >> ltir >> mass >> sfr >> zmet >> tau >> d4000>> qi ;
ssedPhysIn >> age >> luv >> lopt >> lnir >> ltir >> mass >> sfr >> zmet >> tau >> d4000>> qi ;
// read one SED in the binary file
oneSED.readSEDBin(ssedIn);
......@@ -526,7 +523,7 @@ void GalMag::read_SED(){
// emission line fluxes stored into the vector of SED
oneSEDInt.flEm.push_back(val);
// indicate that the emission lines have been computed
oneSEDInt.em='Y';
oneSEDInt.has_emlines = true;
}
// still need to subtract later mag(z=0)
......@@ -584,7 +581,7 @@ Limit the size of the data in memory
*/
GalSED GalMag::reduce_memory(GalSED thisSED){
GalSED reduceSED("intermediate",0,"G");
GalSED reduceSED("intermediate",0);
// Only in the case z>0
if(thisSED.red>1.e-20){
......@@ -657,7 +654,7 @@ void GalMag::comp_mag(){
}
// write the output files
allSED[k].writeLib(outasc,sbinOut,sdatOut,allFlt,magtyp);
allSED[k].writeMag(outasc,sbinOut,sdatOut,allFlt,magtyp);
}
......@@ -683,7 +680,7 @@ void GalMag::print_info(){
cout << endl << "# EM_LINES " << emlines << endl;
cout << "# EM_DISPERSION ";
for (size_t k=0 ; k<fracEm.size() ; k++){cout << fracEm[k] << "," ;};
cout << endl << "# LIB_ASCII " << outasc << endl;
cout << endl << "# LIB_ASCII " << (outasc?"YES":"NO") << endl;
time_t result = time(nullptr);
cout << "# CREATION_DATE " << asctime(std::localtime(&result));
cout << "#############################################" << endl;
......@@ -724,7 +721,7 @@ void QSOMag::print_info(){
for (int k=0 ; k<nextlaw ; k++){cout << modext[k*2] << " "<< modext[k+1] << " " ;};
cout << endl << "# EB_V :" ;
for (int k=0 ; k<nebv ; k++){cout << ebv[k] << " " ;};
cout << "# LIB_ASCII " << outasc << endl;
cout << "# LIB_ASCII " << (outasc?"YES":"NO") << endl;
time_t result = time(nullptr);
cout << "# CREATION_DATE " << asctime(std::localtime(&result));
cout << "#############################################" << endl;
......@@ -778,7 +775,7 @@ void QSOMag::read_SED(){
// Select case which need to be considered (no extinction or extinction in the right model range)
// Remove all cases with extinction not in the right model range
if((ebv[j]<1.e-10 && i==0) || (ebv[j]>0 && oneSED.nummod>=(modext[i*2]+10000) && oneSED.nummod<=(modext[i*2+1]+10000))){
if((ebv[j]<1.e-10 && i==0) || (ebv[j]>0 && oneSED.nummod>=(modext[i*2]) && oneSED.nummod<=(modext[i*2+1]))){
// Generate intermediate Continuum SED, since original one must not change
......@@ -884,7 +881,7 @@ void QSOMag::comp_mag(){
}
// write the output files
allSED[k].writeLib(outasc,sbinOut,sdatOut,allFlt,magtyp);
allSED[k].writeMag(outasc,sbinOut,sdatOut,allFlt,magtyp);
}
......@@ -912,7 +909,7 @@ void StarMag::print_info(){
cout << "# STAR_LIB_IN :" << lepharework + "/lib_bin/" + lib + "(.doc & .bin)" << endl;
cout << "# STAR_LIB_OUT :" << lepharework + "/lib_mag/" + colib + "(.doc & .bin)" << endl;
cout << "# LIB_ASCII " << outasc << endl;
cout << "# LIB_ASCII " << (outasc?"YES":"NO") << endl;
time_t result = time(nullptr);
cout << "# CREATION_DATE " << asctime(std::localtime(&result));
cout << "#############################################" << endl;
......@@ -982,7 +979,7 @@ void StarMag::comp_mag(){
}
// write the output files
(allSED[k]).writeLib(outasc,sbinOut,sdatOut,allFlt,magtyp);
(allSED[k]).writeMag(outasc,sbinOut,sdatOut,allFlt,magtyp);
}
}
......@@ -22,7 +22,8 @@ protected:
string config, typ, verbose;
cosmo lcdm;
string filtFile, magtyp, outasc;
string filtFile, magtyp;
bool outasc;
vector<string> extlaw;
int nextlaw;
vector<double> ebv;
......
......@@ -982,7 +982,7 @@ void onesource::interp_lib(vector<SED*> fulllib,const int imagm, cosmo lcdm){
// Check that the interpolation can be done
// same model in the library around zs
if(SEDa.rec == SEDb.rec && SEDa.ebv == SEDb.ebv && deca!=0){
if(SEDa.is_same_model(SEDb) && deca!=0){
// Linear factor for interpolation
double a=abs((consiz-SEDa.red)/(SEDa.red - SEDb.red));
......@@ -1016,7 +1016,7 @@ void onesource::interp_lib(vector<SED*> fulllib,const int imagm, cosmo lcdm){
//Retrieve the rest-frame magnitudes
int index0=indmin[0];
// go down in the index as long as the redshift >0 and that the model is the same
while((*(fulllib[index0])).red >1.e-20 && (SEDa.rec == (*(fulllib[index0])).rec)){index0--;}
while( (*(fulllib[index0])).red >1.e-20 && SEDa.is_same_model(*fulllib[index0]) ){ index0--; }
if((*(fulllib[index0])).red>1.e-20) cout << " do not find redshift 0 for rest-frame colors. Problem." << endl;
// Store the predicted magnitudes at z=0
for (int k=0; k<imagm; k++){magm0.push_back((*(fulllib[index0])).mag_ori[k]);}
......@@ -1233,7 +1233,7 @@ void onesource::computeEmFlux(vector<SED*> fulllib, cosmo lcdm,vector<opa> opaAl
if(indmin[0]>0){
// Do it only if emission have been generated
if(((fulllib[indmin[0]])->em)[0]!= 'N'){
if(fulllib[indmin[0]]->has_emlines){
// Generate a SED based on the index at z=0: extract the rest-frame SED
ind0=(fulllib[indmin[0]])->index_z0;
......
......@@ -458,21 +458,21 @@ void read_lib(vector<SED*> & libFull, int & ind, int nummodpre[3], const string
SED *oneSED;
if(valc[0] == 'G' || valc[0] == 'g'){
oneSED = new GalSED("bid.dat",0,"G");
oneSED = new GalSED("bid.dat",0);
// If emission lines used, only for galaxies
if(emlines[0]=='P' ||emlines[0]=='E') oneSED->em="YES";
if(emlines[0]=='P' ||emlines[0]=='E') oneSED->has_emlines = true;
indBib=0;
}else if(valc[0] == 'Q' || valc[0] == 'q'){
oneSED = new QSOSED("bid.dat",0,"Q");
oneSED = new QSOSED("bid.dat",0);
indBib=1;
}else if(valc[0] == 'S' || valc[0] == 's'){
oneSED = new StarSED("bid.dat",0,"S");
oneSED = new StarSED("bid.dat",0);
indBib=2;
}else{
throw invalid_argument("There is no such SED type defined: " + valc);
}
// read each SED in the library binary file
oneSED->readLib(slibIn);
oneSED->readMagBin(slibIn);
// fill the vector including the modeled fluxes
oneSED->mag2flux();
// Case in which we add the emission lines (only for galaxies), if option turned on
......
......@@ -5,7 +5,7 @@
import numpy as np
import matplotlib
matplotlib.use('Agg')
matplotlib.use('TkAgg')
import matplotlib.pyplot as plt
from matplotlib.backends.backend_pdf import PdfPages
import sys, getopt
......@@ -29,6 +29,7 @@ for opt, arg in options:
sys.exit()
odev='.'+arg
print('output device will be in '+odev+' format.') #can be png,pdf (if none print on screen)
matplotlib.use('TkAgg')
if opt in ('-c','--context'):
print('context for models: ', arg)
ctx=arg
......
#!/usr/bin/env bash
SCRIPT_DIR="$( cd -- "$( dirname -- "${BASH_SOURCE[0]:-$0}"; )" &> /dev/null && pwd 2> /dev/null; )";
export LEPHAREDIR=$SCRIPT_DIR/..
export LEPHAREWORK=
if [ -z "$1" ]
then
LEPHAREWORK=$SCRIPT_DIR/WORK
else
LEPHAREWORK=$1
fi
export LEPHAREWORK
echo $LEPHAREWORK
mkdir -p $LEPHAREWORK/filt $LEPHAREWORK/lib_bin $LEPHAREWORK/lib_mag $LEPHAREWORK/zphota
export OMP_NUM_THREADS='30'
$LEPHAREDIR/bin/filter -c $LEPHAREDIR/test/COSMOS.para
$LEPHAREDIR/bin/sedtolib -c $LEPHAREDIR/test/COSMOS.para -t S --STAR_SED $LEPHAREDIR/test/STAR_MOD_ALL.list --LIB_ASCII YES
$LEPHAREDIR/bin/mag_gal -c $LEPHAREDIR/test/COSMOS.para -t S --LIB_ASCII YES --STAR_LIB_OUT ALLSTAR_COSMOS
$LEPHAREDIR/bin/sedtolib -c $LEPHAREDIR/test/COSMOS.para -t Q --QSO_SED $LEPHAREDIR/sed/QSO/SALVATO09/AGN_MOD.list
$LEPHAREDIR/bin/mag_gal -c $LEPHAREDIR/test/COSMOS.para -t Q --MOD_EXTINC 0,1000 --EB_V 0.,0.1,0.2,0.3 --EXTINC_LAW SB_calzetti.dat --LIB_ASCII NO --Z_STEP 0.04,0,6 --LIB_ASCII YES
$LEPHAREDIR/bin/sedtolib -c $LEPHAREDIR/test/COSMOS.para -t G --GAL_SED $LEPHAREDIR/test/COSMOS_MOD.list --GAL_LIB LIB_VISTA
$LEPHAREDIR/bin/mag_gal -c $LEPHAREDIR/test/COSMOS.para -t G --GAL_LIB_IN LIB_VISTA --GAL_LIB_OUT VISTA_COSMOS_FREE --MOD_EXTINC 18,26,26,33,26,33,26,33 --EXTINC_LAW SMC_prevot.dat,SB_calzetti.dat,SB_calzetti_bump1.dat,SB_calzetti_bump2.dat --EM_LINES EMP_UV --EM_DISPERSION 0.5,0.75,1.,1.5,2. --Z_STEP 0.04,0,6 --LIB_ASCII YES
cat_out=zphot_short.out
$LEPHAREDIR/bin/zphota -c $LEPHAREDIR/test/COSMOS.para --CAT_IN $LEPHAREDIR/test/COSMOS.in --CAT_OUT $cat_out --ZPHOTLIB VISTA_COSMOS_FREE,ALLSTAR_COSMOS,QSO_COSMOS --ADD_EMLINES 0,100 --AUTO_ADAPT YES --Z_STEP 0.04,0,6 --CAT_LINES 1,100 --SPEC_OUT YES --PARA_OUT $LEPHAREDIR/test/output.para --VERBOSE NO --ZFIX YES
mv $cat_out Id*.spec $LEPHAREWORK/zphota/
#!/usr/bin/env bash
SCRIPT_DIR="$( cd -- "$( dirname -- "${BASH_SOURCE[0]:-$0}"; )" &> /dev/null && pwd 2> /dev/null; )";
export LEPHAREDIR=$SCRIPT_DIR/..
export LEPHAREWORK=
if [ -z "$1" ]
then
LEPHAREWORK=$SCRIPT_DIR/WORK
else
LEPHAREWORK=$1
fi
export LEPHAREWORK
echo $LEPHAREWORK
mkdir -p $LEPHAREWORK/filt $LEPHAREWORK/lib_bin $LEPHAREWORK/lib_mag
export OMP_NUM_THREADS='30'
python $LEPHAREDIR/lephare/filter.py -c $LEPHAREDIR/test/COSMOS.para
python $LEPHAREDIR/lephare/sedtolib.py -c $LEPHAREDIR/test/COSMOS.para --typ S --STAR_SED $LEPHAREDIR/test/STAR_MOD_ALL.list
python $LEPHAREDIR/lephare/mag_gal.py -c $LEPHAREDIR/test/COSMOS.para --typ S --LIB_ASCII YES --STAR_LIB_OUT ALLSTAR_COSMOS
python $LEPHAREDIR/lephare/sedtolib.py -c $LEPHAREDIR/test/COSMOS.para --typ Q --QSO_SED $LEPHAREDIR/sed/QSO/SALVATO09/AGN_MOD.list
python $LEPHAREDIR/lephare/mag_gal.py -c $LEPHAREDIR/test/COSMOS.para --typ Q --MOD_EXTINC 0,1000 --EB_V 0.,0.1,0.2,0.3 --EXTINC_LAW SB_calzetti.dat --LIB_ASCII NO --Z_STEP 0.04,0,6
python $LEPHAREDIR/lephare/sedtolib.py -c $LEPHAREDIR/test/COSMOS.para --typ G --GAL_SED $LEPHAREDIR/test/COSMOS_MOD.list --GAL_LIB LIB_VISTA
python $LEPHAREDIR/lephare/mag_gal.py -c $LEPHAREDIR/test/COSMOS.para --typ G --GAL_LIB_IN LIB_VISTA --GAL_LIB_OUT VISTA_COSMOS_FREE --MOD_EXTINC 18,26,26,33,26,33,26,33 --EXTINC_LAW SMC_prevot.dat,SB_calzetti.dat,SB_calzetti_bump1.dat,SB_calzetti_bump2.dat --EM_LINES EMP_UV --EM_DISPERSION 0.5,0.75,1.,1.5,2. --Z_STEP 0.04,0,6
python $LEPHAREDIR/lephare/zphota.py -c $LEPHAREDIR/test/COSMOS.para --CAT_IN $LEPHAREDIR/test/COSMOS.in --CAT_OUT zphot_short.out --ZPHOTLIB VISTA_COSMOS_FREE,ALLSTAR_COSMOS,QSO_COSMOS --ADD_EMLINES 0,100 --AUTO_ADAPT YES --Z_STEP 0.04,0,6 --CAT_LINES 1,100 --SPEC_OUT YES --PARA_OUT $LEPHAREDIR/test/output.para
  • I am fine with all the modifications here.

    I tried to understand the reason of jtype, nrec, and the convention of nummod (-1* for the stars, +10000 for the QSO). In fact, all of this was to follow the fortran version of the code. But that's true that we don't need this information anymore in the c++ version. One consequence is that the format of some intermediate catalogues changes compared to fortran, but it's not important because it's mainly in the binary files (except the output file created with LIB_ASCII YES, but that's a minor change, easy to spot for the user). So, I think that's fine to simplify the code like this.

    About the use of auto& to read and write the files: that's really efficient. I didn't know about this command.

    I started to look at the python part. Could you change the name of the directory "lephare" in something more meaningful? It was not clear that this directory contains the python codes. Maybe pylephare?

Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment