Manticore  Version 1.5.3
Physics of Molecular Clouds
Namespaces | Classes | Typedefs | Functions | Variables
manticore Namespace Reference

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::stringsetInstruments (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...
 

Detailed Description

Package namespace.

Typedef Documentation

◆ mapDataType

using manticore::mapDataType = typedef std::valarray<float>

Output FITS maps data type (cf. mapDataCode).

Definition at line 12 of file manticore.h.

◆ tableType

using manticore::tableType = typedef std::vector<std::pair<double,double> >

Basic table data.

Definition at line 21 of file manticore.h.

Function Documentation

◆ asctime()

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().

◆ calcStats()

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.

Parameters
[out]statsStatistics storage.
[in]x0Data samples.
[in]nNumber of data samples.

Definition at line 323 of file solve.cc.

References std::sort().

Referenced by findParameters().

◆ ccfCheck()

void manticore::ccfCheck ( )

◆ findMedianModel()

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.

Parameters
[in]pstatsParameters statistics from calcStats().
[in]nParamNumber of parameters.
[in]pdataModel parameter storage array.
[in]nModelsNumber of models (cannot exceed nCols).
[in]nColsNumber of columns per row in pdata.
Returns
Index of optimal model.

Definition at line 354 of file solve.cc.

Referenced by findParameters().

◆ findParameters()

template<unsigned NP>
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).

Template Parameters
NPNumber of (active) parameters.
Parameters
[in]bandDataInput band data.
[in]bandErrorInput band (1-sigma) uncertainties.
[in]tempMapSingle-component temperature map (can be empty).
[in]chi2MapSingle-component reduced \(\chi^2\) map (can be empty).
[in]cliCommand line options and arguments.
[in]gdata0Model parameters.
[in]idThread ID number (0-based).
[in]nThreadsNumber of threads.
Note
When using the –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().

◆ getBand()

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().

◆ getError() [1/2]

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.

Parameters
[in]dataInput image map data.
[in]hduInput image map FITS header.
[in]cliCommand line options and arguments.
[in]llCutout lower-left vertex.
[in]urCutout upper-right vertex.
[in]stCutout stride (cutout disabled if st[0] = 0).
Returns
Image 1-sigma error map.

Definition at line 33 of file process.cc.

References getBand(), mutils::CommandLine::getDouble(), mutils::CommandLine::getOption(), MU_EXCEPTION_IF, and MU_INFO.

◆ getError() [2/2]

mapDataType manticore::getError ( const mapDataType data,
CCfits::PHDU &  hdu,
const mutils::CommandLine cli 
)

Referenced by process().

◆ getSingleMap()

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.

Parameters
[in]cliCommand line options and arguments.
[in]llCutout lower-left vertex.
[in]urCutout upper-right vertex.
[in]stCutout stride (cutout disabled if st[0] = 0).
[in]hduNameMap HDU name to extract.
Returns
Single-temperature reduced \(\chi^2\) map.

Definition at line 953 of file solve.cc.

References std::string::front(), mutils::CommandLine::getOption(), and MU_WARN.

Referenced by solve().

◆ gridSample()

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.

Parameters
[in]dataBand data sample.
[in]errSample uncertainty.
[in]iBandBand index.
[in]iCycleCycle index.
[in]nSampTotal samples per band per cycle.

Definition at line 255 of file solve.cc.

Referenced by findParameters(), and solve().

◆ initSpline()

gsl_spline * manticore::initSpline ( const tableType table,
bool  ln,
unsigned  stride 
)

Initializes cubic spline interpolator.

Note
The input table is not required for interpolation following initialization.
Parameters
[in]tableInterpolation table.
[in]lnWhether to take natural logarithm of the table ordinates.
[in]strideTable data stride (for reducing table length).
Returns
Initialized GSL spline instance.

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().

◆ makeModelMap()

void manticore::makeModelMap ( std::vector< mapDataType > &  bandData,
std::vector< mapDataType > &  bandError,
const GrayData gdata,
unsigned  nParam 
)

Replace input observations with theoretical model.

Parameters
[in,out]bandDataBand data storage.
[in,out]bandErrorBand error storage.
[in]gdataModel parameters.
[in]nParamNumber 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.

Precondition
The input maps must be at least modelMapLen on a side.

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().

◆ printVersion()

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().

◆ process()

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:

  • reading the input FITS images and headers
  • reading or generating input map uncertainties
  • performing input compatibility checks
  • populating the output FITS header and data structures
  • writing the output FITS file to disk
Parameters
[in]cliCommand 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().

◆ resetRefPix()

void manticore::resetRefPix ( CCfits::HDU &  hdu,
const std::vector< long > &  ll 
)

Reset cutout header coordinates.

Parameters
[in,out]hduInput image map FITS header.
[in]llLower-left cutout vertex (cutout disabled if zero).

Definition at line 82 of file process.cc.

Referenced by process(), and setKeys().

◆ robustFit()

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.

Parameters
[in,out]probCeres problem.
[in,out]sumCeres summary.
[in,out]optsCeres options.
[in]costCeres cost function.
[in]pParameters.
[in]dpParameter uncertainties (for refit use; can be null).
[in]pminParameter lower bounds.
[in]pmaxParameter upper bounds.
[in]NPNumber of parameters.
[in]minCostCost below which refit option is ignored.
[in]maxCostCost above which a refit is forced.
[in]refitWhether to force refit.
Returns
Total fitter iterations consumed.

Definition at line 433 of file solve.cc.

References std::max(), std::min(), and trialFit().

Referenced by findParameters().

◆ setInstruments()

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.

Parameters
[in]bandDetector bands.

Definition at line 497 of file solve.cc.

References std::vector::push_back().

Referenced by findParameters().

◆ setKeys()

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.

Parameters
[in,out]hduOutput FITS header.
[in]hdu0Fiducial input image header.
[in]cliCommand line options and arguments.
[in]llLower-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().

◆ solve()

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:

  • T (gas temperature)
  • dT (1-sigma uncertainty in T)
  • N(H2) (N(H2) gas column density)
  • dN(H2) (1-sigma uncertainty in N)
  • Chi2 (fit chi-squared per degree of freedom)
Parameters
[out]outMapOutput result data.
[in,out]outHDUOutput FITS headers.
[in,out]outFITSOutput FITS file.
[in]bandDataInput band data.
[in]bandErrorInput band (1-sigma) uncertainties.
[in]bandHDUInput band FITS headers.
[in]axisImages dimensions.
[in]cliCommand line options and arguments.
[in]llCutout lower-left vertex.
[in]urCutout upper-right vertex.
[in]stCutout 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().

◆ sysSample()

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.

Parameters
[in]dataBand data sample.
[in]iBandBand index.
[in]iCycleCycle index.
[in]instrUnique instruments.
[in]bandsDetector bands.

Definition at line 286 of file solve.cc.

References MU_DEBUG, and std::vector::size().

Referenced by findParameters().

◆ trialFit()

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.

Parameters
[in,out]probCeres problem.
[in,out]sumCeres summary.
[in,out]optsCeres options.
[in]costCeres cost function.
[in]pParameters.
[in]NPNumber of parameters.
[in]p0Previous parameters (can be null).
[in]sum0Previous Ceres summary (can be null).
Returns
Number of fitter iterations.

Definition at line 387 of file solve.cc.

Referenced by robustFit().

Variable Documentation

◆ c_light

constexpr double manticore::c_light = 2.99792458e10

◆ h_Planck

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().

◆ imu_g

const double manticore::imu_g = 1.0/mu_g

Inverse mean molecular weight.

Definition at line 40 of file solve.cc.

Referenced by findParameters().

◆ k_Boltzman

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().

◆ m_H

constexpr double manticore::m_H = 1.6737237e-24

Mass of hygdrogen atom (g).

Definition at line 33 of file manticore.h.

◆ mapDataCode

constexpr int manticore::mapDataCode = FLOAT_IMG

Output FITS maps data code (cf. mapDataType).

Definition at line 15 of file manticore.h.

Referenced by process(), and solve().

◆ modelMapLen

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().

◆ mu_g

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().

◆ reldate

constexpr const char* const manticore::reldate = "Jan 15 2020"

Release date (if any).

Definition at line 7 of file version.h.

Referenced by printVersion().

◆ version

constexpr const char* const manticore::version = "1.5.3"

Version string.

Definition at line 6 of file version.h.

Referenced by printVersion(), and setKeys().