Manticore
Version 1.5.3
Physics of Molecular Clouds
|
Package namespace. More...
Namespaces | |
detector | |
Detector support. | |
dust | |
Dust models. | |
graybody | |
Detector support. | |
Classes | |
class | Detector |
Detector base class. More... | |
class | Dust |
Dust model. More... | |
class | Fcost_gray |
Ceres solver graybody cost model. More... | |
class | Graybody |
Graybody emission model. More... | |
struct | GrayData |
Observation fitting data structure. More... | |
class | PACS |
PACS detector class. More... | |
class | SPIRE |
SPIRE detector class. More... | |
Typedefs | |
using | mapDataType = std::valarray< float > |
Output FITS maps data type (cf. mapDataCode). More... | |
using | tableType = std::vector< std::pair< double, double > > |
Basic table data. More... | |
Functions | |
std::string | asctime () |
Current time a la std::asctime(). More... | |
void | calcStats (double stats[], const double x0[], unsigned n) |
Computes median (and arithmetic) statistics on an input vector. More... | |
void | ccfCheck () |
unsigned | findMedianModel (const double pstats[][6], unsigned nParam, const double *pdata, unsigned nModels, unsigned nCols) |
Finds optimal model from a set of alternatives. More... | |
template<unsigned NP> | |
void | findParameters (const std::vector< mapDataType > &bandData, const std::vector< mapDataType > &bandError, const mapDataType &tempMap, const mapDataType &chi2Map, const mu::CommandLine &cli, const GrayData &gdata0, unsigned id, unsigned nThreads) |
Finds best-fit model parameters. More... | |
int | getBand (const CCfits::HDU &hdu) |
Returns nominal band center (microns) from an HDU. More... | |
mapDataType | getError (const mapDataType &data, CCfits::PHDU &hdu, const mu::CommandLine &cli, const std::vector< long > &ll, const std::vector< long > &ur, const std::vector< long > &st) |
Generates image error map. More... | |
mapDataType | getError (const mapDataType &data, CCfits::PHDU &hdu, const mutils::CommandLine &cli) |
mapDataType | getSingleMap (const mu::CommandLine &cli, const std::vector< long > &ll, const std::vector< long > &ur, const std::vector< long > &st, const std::string &hduName) |
Read a map plane from the --single-temp FITS file. More... | |
double | gridSample (double data, double err, unsigned iBand, unsigned iCycle, unsigned nSamp) |
Sample random data errors on a fixed grid. More... | |
gsl_spline * | initSpline (const tableType &table, bool ln, unsigned stride) |
Initializes cubic spline interpolator. More... | |
void | makeModelMap (std::vector< mapDataType > &bandData, std::vector< mapDataType > &bandError, const GrayData &gdata, unsigned nParam) |
Replace input observations with theoretical model. More... | |
void | printVersion () |
Prints version information to the log. More... | |
void | process (const mu::CommandLine &cli) |
Processes multi-wavelength maps into temperature/column density maps. More... | |
void | resetRefPix (CCfits::HDU &hdu, const std::vector< long > &ll) |
Reset cutout header coordinates. More... | |
int | robustFit (ceres::Problem &prob, ceres::Solver::Summary &sum, ceres::Solver::Options &opts, const ceres::CostFunction *cost, double p[], const double dp[], const double pmin[], const double pmax[], unsigned NP, double minCost=0.25, double maxCost=25.0, bool refit=false) |
Non-linear fit with robustness checking. More... | |
std::vector< std::string > | setInstruments (const std::vector< const Detector * > &band) |
Set unique instrument types. More... | |
void | setKeys (CCfits::ExtHDU &hdu, const CCfits::PHDU &hdu0, const mu::CommandLine &cli, const std::vector< long > &ll) |
Sets header keywords for an output map. More... | |
void | solve (std::vector< mapDataType > &outMap, std::vector< CCfits::ExtHDU * > &outHDU, CCfits::FITS *outFITS, const std::vector< mapDataType > &bandData, const std::vector< mapDataType > &bandError, const std::vector< CCfits::PHDU * > &bandHDU, const std::vector< long > &axis, const mu::CommandLine &cli, const std::vector< long > &ll, const std::vector< long > &ur, const std::vector< long > &st) |
Performs least-squares fitting of band SEDs. More... | |
double | sysSample (double data, unsigned iBand, unsigned iCycle, const std::vector< std::string > &instr, const std::vector< const Detector * > &bands) |
Sample systematic data errors on a fixed grid. More... | |
int | trialFit (ceres::Problem &prob, ceres::Solver::Summary &sum, ceres::Solver::Options &opts, const ceres::CostFunction *cost, double p[], unsigned NP, const double p0[]=nullptr, const ceres::Solver::Summary *sum0=nullptr) |
Trial fit with rejection check. More... | |
Variables | |
constexpr double | c_light = 2.99792458e10 |
Speed of light (cm/s). More... | |
constexpr double | h_Planck = 6.626070040e-27 |
Planck's constant (erg*s). More... | |
const double | imu_g = 1.0/mu_g |
Inverse mean molecular weight. More... | |
constexpr double | k_Boltzman = 1.38064852e-16 |
Planck's constant (erg/K). More... | |
constexpr double | m_H = 1.6737237e-24 |
Mass of hygdrogen atom (g). More... | |
constexpr int | mapDataCode = FLOAT_IMG |
Output FITS maps data code (cf. mapDataType). More... | |
constexpr long | modelMapLen = 961L |
One-dimensional length of model map image. More... | |
const double | mu_g = (2.75*m_H) |
Mean molecular weight. More... | |
constexpr const char *const | reldate = "Jan 15 2020" |
Release date (if any). More... | |
constexpr const char *const | version = "1.5.3" |
Version string. More... | |
Package namespace.
using manticore::mapDataType = typedef std::valarray<float> |
Output FITS maps data type (cf. mapDataCode).
Definition at line 12 of file manticore.h.
using manticore::tableType = typedef std::vector<std::pair<double,double> > |
Basic table data.
Definition at line 21 of file manticore.h.
std::string manticore::asctime | ( | ) |
Current time a la std::asctime().
Definition at line 182 of file solve.cc.
References std::string::pop_back(), and std::time().
Referenced by findParameters().
void manticore::calcStats | ( | double | stats[], |
const double | x0[], | ||
unsigned | n | ||
) |
Computes median (and arithmetic) statistics on an input vector.
The input vector is copied, sorted and analyzed to determine the median, quantile-based deviation, arithmtic mean, sample standard deviation, minimum and maximum, which are written to stats in that order. Hence stats must contain at least six elements.
[out] | stats | Statistics storage. |
[in] | x0 | Data samples. |
[in] | n | Number of data samples. |
Definition at line 323 of file solve.cc.
References std::sort().
Referenced by findParameters().
void manticore::ccfCheck | ( | ) |
Definition at line 17 of file manticore.cc.
References c_light, manticore::Detector::ccf(), manticore::Graybody::dF_dS(), manticore::Graybody::dF_dT(), manticore::Graybody::F(), manticore::Detector::freqWidth(), manticore::PACS::init(), manticore::Graybody::Inu(), manticore::Dust::setModel(), and mutils::Log::stream.
Referenced by main().
unsigned manticore::findMedianModel | ( | const double | pstats[][6], |
unsigned | nParam, | ||
const double * | pdata, | ||
unsigned | nModels, | ||
unsigned | nCols | ||
) |
Finds optimal model from a set of alternatives.
The optimal model is considered the one with minimum deviation from the median model parameters, taking also its \(\chi^2\) into account.
[in] | pstats | Parameters statistics from calcStats(). |
[in] | nParam | Number of parameters. |
[in] | pdata | Model parameter storage array. |
[in] | nModels | Number of models (cannot exceed nCols). |
[in] | nCols | Number of columns per row in pdata. |
Definition at line 354 of file solve.cc.
Referenced by findParameters().
void manticore::findParameters | ( | const std::vector< mapDataType > & | bandData, |
const std::vector< mapDataType > & | bandError, | ||
const mapDataType & | tempMap, | ||
const mapDataType & | chi2Map, | ||
const mu::CommandLine & | cli, | ||
const GrayData & | gdata0, | ||
unsigned | id, | ||
unsigned | nThreads | ||
) |
Finds best-fit model parameters.
Solves the non-linear least-squares problem to find the best-fit temperature and surface density for up to two gas components. The one-component model consists of a single "hot"/"halo" slab of gas. The two-component model consists of a "cold"/"core" slab sandwiched between two identical "hot"/"halo" slabs, for a total of four parameters (two temperatures and surface densities). Optionally, the halo temperature can be prescribed and held fixed (NP = 1 or NP = 3).
NP | Number of (active) parameters. |
[in] | bandData | Input band data. |
[in] | bandError | Input band (1-sigma) uncertainties. |
[in] | tempMap | Single-component temperature map (can be empty). |
[in] | chi2Map | Single-component reduced \(\chi^2\) map (can be empty). |
[in] | cli | Command line options and arguments. |
[in] | gdata0 | Model parameters. |
[in] | id | Thread ID number (0-based). |
[in] | nThreads | Number of threads. |
–sys-error
option, the random error maps should probably be systematically rescaled the same as the data, but this has not been implemented. Definition at line 540 of file solve.cc.
References asctime(), manticore::GrayData::band, calcStats(), mutils::Log::DEBUG, std::vector::emplace_back(), findMedianModel(), mutils::CommandLine::getDouble(), mutils::CommandLine::getInteger(), mutils::CommandLine::getOption(), manticore::GrayData::gray, gridSample(), imu_g, std::isfinite(), mutils::Log::level, std::max(), std::min(), MU_DEBUG, MU_ERROR, mu_g, MU_INFO, MU_INFO_IF, MU_WARN, manticore::GrayData::nBands, manticore::GrayData::outMap, robustFit(), setInstruments(), std::vector::size(), mutils::Log::stream, and sysSample().
int manticore::getBand | ( | const CCfits::HDU & | hdu | ) |
Returns nominal band center (microns) from an HDU.
Definition at line 105 of file manticore.cc.
Referenced by getError(), process(), and solve().
mapDataType manticore::getError | ( | const mapDataType & | data, |
CCfits::PHDU & | hdu, | ||
const mu::CommandLine & | cli, | ||
const std::vector< long > & | ll, | ||
const std::vector< long > & | ur, | ||
const std::vector< long > & | st | ||
) |
Generates image error map.
If any --error-map
options were supplied and an entry matching the BAND and TELESCOP keywords of the input is found, the error map is read from the corresponding image. Otherwise, the value of --error-pct
(default: 3.0) is multiplied by data to create the error map.
[in] | data | Input image map data. |
[in] | hdu | Input image map FITS header. |
[in] | cli | Command line options and arguments. |
[in] | ll | Cutout lower-left vertex. |
[in] | ur | Cutout upper-right vertex. |
[in] | st | Cutout stride (cutout disabled if st[0] = 0). |
Definition at line 33 of file process.cc.
References getBand(), mutils::CommandLine::getDouble(), mutils::CommandLine::getOption(), MU_EXCEPTION_IF, and MU_INFO.
mapDataType manticore::getError | ( | const mapDataType & | data, |
CCfits::PHDU & | hdu, | ||
const mutils::CommandLine & | cli | ||
) |
Referenced by process().
mapDataType manticore::getSingleMap | ( | const mu::CommandLine & | cli, |
const std::vector< long > & | ll, | ||
const std::vector< long > & | ur, | ||
const std::vector< long > & | st, | ||
const std::string & | hduName | ||
) |
Read a map plane from the --single-temp
FITS file.
If the --single-temp
option was given, extract the HDU named hduName from the FITS file.
A 'T' map can be used to define the halo temperature for 1- or 3-parameter fits.
A 'Chi2/DOF' map can be used in masking evaluation of the stage 2 (3- or 4-parameter) fitter for individual pixels. The corresponding cut-off is specified via the --chi2-max
option.
[in] | cli | Command line options and arguments. |
[in] | ll | Cutout lower-left vertex. |
[in] | ur | Cutout upper-right vertex. |
[in] | st | Cutout stride (cutout disabled if st[0] = 0). |
[in] | hduName | Map HDU name to extract. |
Definition at line 953 of file solve.cc.
References std::string::front(), mutils::CommandLine::getOption(), and MU_WARN.
Referenced by solve().
double manticore::gridSample | ( | double | data, |
double | err, | ||
unsigned | iBand, | ||
unsigned | iCycle, | ||
unsigned | nSamp | ||
) |
Sample random data errors on a fixed grid.
To empirically explore the effects of data uncertainty on parameter results, this method returns an effective data value using the sample index to cycle through \( data \pm err \) values for each band. For this, bit iBand of iCycle is used as the sign bit of error offset (nSamp = 2). For nSamp = 3, the unaltered data sample is also included in the cycle.
[in] | data | Band data sample. |
[in] | err | Sample uncertainty. |
[in] | iBand | Band index. |
[in] | iCycle | Cycle index. |
[in] | nSamp | Total samples per band per cycle. |
Definition at line 255 of file solve.cc.
Referenced by findParameters(), and solve().
gsl_spline * manticore::initSpline | ( | const tableType & | table, |
bool | ln, | ||
unsigned | stride | ||
) |
Initializes cubic spline interpolator.
table
is not required for interpolation following initialization.[in] | table | Interpolation table. |
[in] | ln | Whether to take natural logarithm of the table ordinates. |
[in] | stride | Table data stride (for reducing table length). |
Definition at line 81 of file manticore.cc.
References std::vector::push_back(), and std::vector::size().
Referenced by manticore::Detector::init(), and manticore::Dust::setModel().
void manticore::makeModelMap | ( | std::vector< mapDataType > & | bandData, |
std::vector< mapDataType > & | bandError, | ||
const GrayData & | gdata, | ||
unsigned | nParam | ||
) |
Replace input observations with theoretical model.
[in,out] | bandData | Band data storage. |
[in,out] | bandError | Band error storage. |
[in] | gdata | Model parameters. |
[in] | nParam | Number of parameters. |
The model varies the parameters Tc, Sc, Th, Sh on a regular grid, core parameters varying along each row (temperature first) and the halo parameters varying along the columns.
Definition at line 208 of file solve.cc.
References manticore::GrayData::band, manticore::GrayData::dnu0, manticore::GrayData::gray, modelMapLen, MU_DEBUG, mu_g, and std::vector::size().
Referenced by solve().
void manticore::printVersion | ( | ) |
Prints version information to the log.
Definition at line 94 of file manticore.cc.
References std::string::c_str(), reldate, mutils::Log::stream, and version.
Referenced by main().
void manticore::process | ( | const mu::CommandLine & | cli | ) |
Processes multi-wavelength maps into temperature/column density maps.
This routine handles the higher-level management required to generate mass and temperature images from the input continuum data, including:
[in] | cli | Command line options and arguments. |
Definition at line 104 of file process.cc.
References mutils::CommandLine::arguments(), std::vector::at(), std::vector::back(), std::vector::emplace_back(), std::vector::front(), std::unique_ptr::get(), getBand(), getError(), mutils::CommandLine::getOption(), mutils::CommandLine::getString(), std::string::insert(), mapDataCode, manticore::dust::modelMap, modelMapLen, MU_ERROR_IF, MU_EXCEPTION_IF, MU_INFO, mutils::nullStr, std::vector::push_back(), resetRefPix(), setKeys(), std::vector::size(), solve(), and std::to_string().
Referenced by main().
void manticore::resetRefPix | ( | CCfits::HDU & | hdu, |
const std::vector< long > & | ll | ||
) |
Reset cutout header coordinates.
[in,out] | hdu | Input image map FITS header. |
[in] | ll | Lower-left cutout vertex (cutout disabled if zero). |
Definition at line 82 of file process.cc.
int manticore::robustFit | ( | ceres::Problem & | prob, |
ceres::Solver::Summary & | sum, | ||
ceres::Solver::Options & | opts, | ||
const ceres::CostFunction * | cost, | ||
double | p[], | ||
const double | dp[], | ||
const double | pmin[], | ||
const double | pmax[], | ||
unsigned | NP, | ||
double | minCost = 0.25 , |
||
double | maxCost = 25.0 , |
||
bool | refit = false |
||
) |
Non-linear fit with robustness checking.
If the initial fit fails by one of several criteria, or if explicitly requested, the parameters are perturbed and refit; if more successful, the refit replaces the original solution. Parameter perturbations are set by the dp[] uncertainties, with signs interleaved between component parameters to account for typical covariance between them. If dp is null, a relative error of 10% is used.
[in,out] | prob | Ceres problem. |
[in,out] | sum | Ceres summary. |
[in,out] | opts | Ceres options. |
[in] | cost | Ceres cost function. |
[in] | p | Parameters. |
[in] | dp | Parameter uncertainties (for refit use; can be null). |
[in] | pmin | Parameter lower bounds. |
[in] | pmax | Parameter upper bounds. |
[in] | NP | Number of parameters. |
[in] | minCost | Cost below which refit option is ignored. |
[in] | maxCost | Cost above which a refit is forced. |
[in] | refit | Whether to force refit. |
Definition at line 433 of file solve.cc.
References std::max(), std::min(), and trialFit().
Referenced by findParameters().
std::vector<std::string> manticore::setInstruments | ( | const std::vector< const Detector * > & | band | ) |
Set unique instrument types.
This is used to associate correlated systematic errors with the appropriate set of detector bands. It is assumed band names follow the INSTRUMENT-BAND
convention, e.g., PACS-160
.
[in] | band | Detector bands. |
Definition at line 497 of file solve.cc.
References std::vector::push_back().
Referenced by findParameters().
void manticore::setKeys | ( | CCfits::ExtHDU & | hdu, |
const CCfits::PHDU & | hdu0, | ||
const mu::CommandLine & | cli, | ||
const std::vector< long > & | ll | ||
) |
Sets header keywords for an output map.
Combines information from the input FITS headers and the command line arguments to populate keywords in the output FITS file.
[in,out] | hdu | Output FITS header. |
[in] | hdu0 | Fiducial input image header. |
[in] | cli | Command line options and arguments. |
[in] | ll | Lower-left cutout vertex (ignored if zero). |
Definition at line 241 of file process.cc.
References std::string::find(), mutils::CommandLine::getOption(), MU_EXCEPTION, resetRefPix(), and version.
Referenced by process().
void manticore::solve | ( | std::vector< mapDataType > & | outMap, |
std::vector< CCfits::ExtHDU * > & | outHDU, | ||
CCfits::FITS * | outFITS, | ||
const std::vector< mapDataType > & | bandData, | ||
const std::vector< mapDataType > & | bandError, | ||
const std::vector< CCfits::PHDU * > & | bandHDU, | ||
const std::vector< long > & | axis, | ||
const mu::CommandLine & | cli, | ||
const std::vector< long > & | ll, | ||
const std::vector< long > & | ur, | ||
const std::vector< long > & | st | ||
) |
Performs least-squares fitting of band SEDs.
Given input band continuum fluxes and uncertainties, calculates the best-fit temperature and column density maps (including uncertainties and reduced chi-square information).
For one component the output FITS file will contain the following named image extensions:
[out] | outMap | Output result data. |
[in,out] | outHDU | Output FITS headers. |
[in,out] | outFITS | Output FITS file. |
[in] | bandData | Input band data. |
[in] | bandError | Input band (1-sigma) uncertainties. |
[in] | bandHDU | Input band FITS headers. |
[in] | axis | Images dimensions. |
[in] | cli | Command line options and arguments. |
[in] | ll | Cutout lower-left vertex. |
[in] | ur | Cutout upper-right vertex. |
[in] | st | Cutout stride (cutout disabled if st[0] = 0). |
Definition at line 1006 of file solve.cc.
References std::vector::back(), bandData1, std::string::c_str(), std::vector::emplace_back(), std::string::find(), getBand(), mutils::CommandLine::getDouble(), mutils::CommandLine::getInteger(), mutils::CommandLine::getOption(), getSingleMap(), mutils::CommandLine::getString(), gridSample(), makeModelMap(), mapDataCode, std::max(), manticore::dust::modelMap, MU_INFO, MU_INFO_IF, manticore::Dust::name(), mutils::npos, std::vector::push_back(), manticore::Dust::rho(), std::vector::size(), std::stod(), and std::to_string().
Referenced by process().
double manticore::sysSample | ( | double | data, |
unsigned | iBand, | ||
unsigned | iCycle, | ||
const std::vector< std::string > & | instr, | ||
const std::vector< const Detector * > & | bands | ||
) |
Sample systematic data errors on a fixed grid.
Detector objects contain information on instrument systematic errors, both correlated and uncorrelated between bands. This routine cycles through all possibilities of \( data \pm err \)—up to \(2^{NB+NI}\), where NB
is the number of bands and NI
the number of unique instruments between them. If any of the systematic errors is zero, the number of possibilities is reduced. This is indicated by setting the output fluxes to NAN
.
[in] | data | Band data sample. |
[in] | iBand | Band index. |
[in] | iCycle | Cycle index. |
[in] | instr | Unique instruments. |
[in] | bands | Detector bands. |
Definition at line 286 of file solve.cc.
References MU_DEBUG, and std::vector::size().
Referenced by findParameters().
int manticore::trialFit | ( | ceres::Problem & | prob, |
ceres::Solver::Summary & | sum, | ||
ceres::Solver::Options & | opts, | ||
const ceres::CostFunction * | cost, | ||
double | p[], | ||
unsigned | NP, | ||
const double | p0[] = nullptr , |
||
const ceres::Solver::Summary * | sum0 = nullptr |
||
) |
Trial fit with rejection check.
Performs a trial fit and compares it to a previous fit (if provided). Restores the old fit if that was better.
[in,out] | prob | Ceres problem. |
[in,out] | sum | Ceres summary. |
[in,out] | opts | Ceres options. |
[in] | cost | Ceres cost function. |
[in] | p | Parameters. |
[in] | NP | Number of parameters. |
[in] | p0 | Previous parameters (can be null). |
[in] | sum0 | Previous Ceres summary (can be null). |
Definition at line 387 of file solve.cc.
Referenced by robustFit().
constexpr double manticore::c_light = 2.99792458e10 |
Speed of light (cm/s).
Definition at line 24 of file manticore.h.
Referenced by manticore::Detector::calcWidth(), ccfCheck(), manticore::SPIRE::init(), manticore::PACS::init(), manticore::Detector::init(), and main().
constexpr double manticore::h_Planck = 6.626070040e-27 |
Planck's constant (erg*s).
Definition at line 27 of file manticore.h.
Referenced by manticore::Graybody::Bnu(), and manticore::Graybody::dBnu_dT().
const double manticore::imu_g = 1.0/mu_g |
Inverse mean molecular weight.
Definition at line 40 of file solve.cc.
Referenced by findParameters().
constexpr double manticore::k_Boltzman = 1.38064852e-16 |
Planck's constant (erg/K).
Definition at line 30 of file manticore.h.
Referenced by manticore::Graybody::Bnu(), and manticore::Graybody::dBnu_dT().
constexpr double manticore::m_H = 1.6737237e-24 |
Mass of hygdrogen atom (g).
Definition at line 33 of file manticore.h.
constexpr int manticore::mapDataCode = FLOAT_IMG |
Output FITS maps data code (cf. mapDataType).
Definition at line 15 of file manticore.h.
constexpr long manticore::modelMapLen = 961L |
One-dimensional length of model map image.
Definition at line 18 of file manticore.h.
Referenced by makeModelMap(), and process().
const double manticore::mu_g = (2.75*m_H) |
Mean molecular weight.
Definition at line 39 of file solve.cc.
Referenced by findParameters(), and makeModelMap().
constexpr const char* const manticore::reldate = "Jan 15 2020" |
constexpr const char* const manticore::version = "1.5.3" |
Version string.
Definition at line 6 of file version.h.
Referenced by printVersion(), and setKeys().