Manticore  Version 2.0alpha
Physics of Molecular Clouds
manticore.h
Go to the documentation of this file.
1 #ifndef MANTICORE_H
2 #define MANTICORE_H
3 
4 #include <CCfits/CCfits>
5 #include <gsl/gsl_spline.h>
6 #include <mutils/CommandLine.h>
7 
8 namespace manticore {
9 
11 using mapDataType = std::valarray<float>;
12 
14 constexpr int mapDataCode = FLOAT_IMG;
15 
17 constexpr long modelMapLen = 961L;
18 
20 using tableType = std::vector<std::pair<double,double>>;
21 
23 constexpr double c_light = 2.99792458e10;
24 
26 constexpr double h_Planck = 6.626070040e-27;
27 
29 constexpr double k_Boltzman = 1.38064852e-16;
30 
32 constexpr double m_H = 1.6737237e-24;
33 
34 extern mapDataType
35  getError(const mapDataType &data, CCfits::PHDU &hdu,
36  const mutils::CommandLine &cli);
37 
38 extern gsl_spline
39  *initSpline(const tableType &table, bool ln = false);
40 
41 extern void
42  process(const mutils::CommandLine &cli);
43 
44 extern void
45  setKeys(CCfits::ExtHDU &hdu, const CCfits::PHDU &hdu0);
46 
47 extern void
48  solve(std::vector<mapDataType> &outMap,
49  std::vector<CCfits::ExtHDU *> &outHDU,
50  CCfits::FITS *outFITS,
51  const std::vector<mapDataType> &bandData,
52  const std::vector<mapDataType> &bandError,
53  const std::vector<CCfits::PHDU*> &bandHDU,
54  const mutils::CommandLine &cli);
55 
56 extern int
57  getBand(const CCfits::HDU &hdu);
58 
59 } // namespace manticore
60 
61 #endif // MANTICORE_H
void process(const mutils::CommandLine &cli)
Processes multi-wavelength maps into temperature/column density maps.
Definition: process.cc:82
void setKeys(CCfits::ExtHDU &hdu, const CCfits::PHDU &hdu0)
Sets header keywords for an output map.
Definition: process.cc:196
constexpr double h_Planck
Planck&#39;s constant (erg*s).
Definition: manticore.h:26
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 mutils::CommandLine &cli)
Performs least-squares fitting of band SEDs.
Definition: solve.cc:635
gsl_spline * initSpline(const tableType &table, bool ln)
Initializes cubic spline interpolator.
Definition: manticore.cc:74
constexpr double k_Boltzman
Planck&#39;s constant (erg/K).
Definition: manticore.h:29
Command line options and arguments.
Definition: CommandLine.h:47
Package namespace.
Definition: Detector.cc:13
constexpr double m_H
Mass of hygdrogen atom (g).
Definition: manticore.h:32
Declares the CommandLine class.
constexpr double c_light
Speed of light (cm/s).
Definition: manticore.h:23
mapDataType getError(const mapDataType &data, CCfits::PHDU &hdu, const mutils::CommandLine &cli)
Generates image error map.
Definition: process.cc:30
int getBand(const CCfits::HDU &hdu)
Returns nominal band center (microns) from an HDU.
Definition: manticore.cc:97
constexpr long modelMapLen
One-dimensional length of model map image.
Definition: manticore.h:17
std::vector< std::pair< double, double > > tableType
Basic table data.
Definition: manticore.h:20
constexpr int mapDataCode
Output FITS maps data code (cf. mapDataType).
Definition: manticore.h:14
std::valarray< float > mapDataType
Output FITS maps data type (cf. mapDataCode).
Definition: manticore.h:11