43 using namespace catalogue;
44 using namespace chainmesh;
46 using namespace pairs;
47 using namespace measure;
48 using namespace twopt;
56 if (!compute_extra_info)
57 m_dd = (binType==BinType::_logarithmic_) ? move(Pair::Create(PairType::_angular_log_, PairInfo::_standard_, thetaMin, thetaMax, nbins, shift, angularUnits, angularWeight))
58 : move(Pair::Create(PairType::_angular_lin_, PairInfo::_standard_, thetaMin, thetaMax, nbins, shift, angularUnits, angularWeight));
60 m_dd = (binType==BinType::_logarithmic_) ? move(Pair::Create(PairType::_angular_log_, PairInfo::_extra_, thetaMin, thetaMax, nbins, shift, angularUnits, angularWeight))
61 : move(Pair::Create(PairType::_angular_lin_, PairInfo::_extra_, thetaMin, thetaMax, nbins, shift, angularUnits, angularWeight));
63 m_rr = (binType==BinType::_logarithmic_) ? move(Pair::Create(PairType::_angular_log_, PairInfo::_standard_, thetaMin, thetaMax, nbins, shift, angularUnits))
64 : move(Pair::Create(PairType::_angular_lin_, PairInfo::_extra_, thetaMin, thetaMax, nbins, shift, angularUnits));
66 m_dr = (binType==BinType::_logarithmic_) ? move(Pair::Create(PairType::_angular_log_, PairInfo::_standard_, thetaMin, thetaMax, nbins, shift, angularUnits))
67 : move(Pair::Create(PairType::_angular_lin_, PairInfo::_extra_, thetaMin, thetaMax, nbins, shift, angularUnits));
76 if (!compute_extra_info)
77 m_dd = (binType==BinType::_logarithmic_) ? move(Pair::Create(PairType::_angular_log_, PairInfo::_standard_, thetaMin, thetaMax, binSize, shift, angularUnits, angularWeight))
78 : move(Pair::Create(PairType::_angular_lin_, PairInfo::_standard_, thetaMin, thetaMax, binSize, shift, angularUnits, angularWeight));
80 m_dd = (binType==BinType::_logarithmic_) ? move(Pair::Create(PairType::_angular_log_, PairInfo::_extra_, thetaMin, thetaMax, binSize, shift, angularUnits, angularWeight))
81 : move(Pair::Create(PairType::_angular_lin_, PairInfo::_extra_, thetaMin, thetaMax, binSize, shift, angularUnits, angularWeight));
83 m_rr = (binType==BinType::_logarithmic_) ? move(Pair::Create(PairType::_angular_log_, PairInfo::_standard_, thetaMin, thetaMax, binSize, shift, angularUnits))
84 : move(Pair::Create(PairType::_angular_lin_, PairInfo::_standard_, thetaMin, thetaMax, binSize, shift, angularUnits));
86 m_dr = (binType==BinType::_logarithmic_) ? move(Pair::Create(PairType::_angular_log_, PairInfo::_standard_, thetaMin, thetaMax, binSize, shift, angularUnits))
87 : move(Pair::Create(PairType::_angular_lin_, PairInfo::_standard_, thetaMin, thetaMax, binSize, shift, angularUnits));
96 m_dataset->read(dir+file);
105 vector<double> xx = m_dataset->xx();
107 checkDim(xx, m_dd->nbins(),
"theta");
109 string header =
"[1] angular separation at the bin centre # [2] angular two-point correlation function # [3] error";
110 if (m_compute_extra_info) header +=
" # [4] mean angular separation # [5] standard deviation of the distribution of angular separations # [6] mean redshift # [7] standard deviation of the redshift distribution";
112 m_dataset->write(dir, file, header, 5, rank);
119 void cbl::measure::twopt::TwoPointCorrelation1D_angular::measure (
const ErrorType errorType,
const string dir_output_pairs,
const vector<string> dir_input_pairs,
const string dir_output_resample,
int nMocks,
bool count_dd,
const bool count_rr,
const bool count_dr,
const bool tcount,
const Estimator estimator,
const double fact,
const int seed)
122 case(ErrorType::_Poisson_):
123 measurePoisson(dir_output_pairs, dir_input_pairs, count_dd, count_rr, count_dr, tcount, estimator, fact);
125 case(ErrorType::_Jackknife_):
126 measureJackknife(dir_output_pairs, dir_input_pairs, dir_output_resample, count_dd, count_rr, count_dr, tcount, estimator, fact);
128 case(ErrorType::_Bootstrap_):
129 measureBootstrap(nMocks, dir_output_pairs, dir_input_pairs, dir_output_resample, count_dd, count_rr, count_dr, tcount, estimator, fact, seed);
132 ErrorCBL(
"unknown type of error!",
"measure",
"TwoPointCorrelation1D_angular.cpp");
144 count_allPairs(m_twoPType, dir_output_pairs, dir_input_pairs, count_dd, count_rr, count_dr, tcount, estimator, fact);
149 if (estimator==Estimator::_natural_)
150 m_dataset = correlation_NaturalEstimator(m_dd, m_rr);
151 else if (estimator==Estimator::_LandySzalay_)
152 m_dataset = correlation_LandySzalayEstimator(m_dd, m_rr, m_dr);
154 ErrorCBL(
"the chosen estimator is not implemented!",
"measurePoisson",
"TwoPointCorrelation1D_angular.cpp");
164 string mkdir =
"mkdir -p "+dir_output_resample;
165 if (system(mkdir.c_str())) {}
168 vector<long> region_list = m_data->region_list();
169 size_t nRegions = region_list.size();
171 vector<vector<double> > xi_SubSamples, covariance;
173 vector<shared_ptr<Pair> > dd_regions, rr_regions, dr_regions;
175 count_allPairs_region(dd_regions, rr_regions, dr_regions, m_twoPType, dir_output_pairs, dir_input_pairs, count_dd, count_rr, count_dr, tcount, estimator, fact);
177 vector<shared_ptr<Data> > data_SS = (estimator==Estimator::_natural_) ? XiJackknife(dd_regions, rr_regions) : XiJackknife(dd_regions, rr_regions, dr_regions);
179 for (
size_t i=0; i<nRegions; i++) {
183 string header =
"[1] angular separation at the bin centre # [2] angular two-point correlation function # [3] error";
184 if (m_compute_extra_info) header +=
" # [4] mean separation # [5] standard deviation of the separation distribution # [6] mean redshift # [7] standard deviation of the redshift distribution";
185 data_SS[i]->write(dir_output_resample, file, header, 10, 0);
188 vector<double> dd; data_SS[i]->get_data(dd);
190 xi_SubSamples.push_back(dd);
195 if (estimator==Estimator::_natural_)
196 m_dataset = correlation_NaturalEstimator(m_dd, m_rr);
197 else if (estimator==Estimator::_LandySzalay_)
198 m_dataset = correlation_LandySzalayEstimator(m_dd, m_rr, m_dr);
200 ErrorCBL(
"the chosen estimator is not implemented!",
"measureJackknife",
"TwoPointCorrelation1D_angular.cpp");
202 m_dataset->set_covariance(covariance);
210 void cbl::measure::twopt::TwoPointCorrelation1D_angular::measureBootstrap (
const int nMocks,
const string dir_output_pairs,
const vector<string> dir_input_pairs,
const string dir_output_resample,
const bool count_dd,
const bool count_rr,
const bool count_dr,
const bool tcount,
const Estimator estimator,
const double fact,
const int seed)
213 ErrorCBL(
"the number of mocks must be >0!",
"measureBootstrap",
"TwoPointCorrelation1D_angular.cpp");
216 string mkdir =
"mkdir -p "+dir_output_resample;
217 if (system(mkdir.c_str())) {}
220 vector<vector<double> > xi_SubSamples, covariance;
222 vector<shared_ptr<Pair> > dd_regions, rr_regions, dr_regions;
224 count_allPairs_region(dd_regions, rr_regions, dr_regions, m_twoPType, dir_output_pairs, dir_input_pairs, count_dd, count_rr, count_dr, tcount, estimator, fact);
226 vector<shared_ptr<Data> > data_SS = (estimator==Estimator::_natural_) ? XiBootstrap(nMocks, dd_regions, rr_regions, seed) : XiBootstrap(nMocks, dd_regions, rr_regions, dr_regions, seed);
228 for (
int i=0; i<nMocks; i++) {
232 string header =
"[1] angular separation at the bin centre # [2] angular two-point correlation function # [3] error";
233 if (m_compute_extra_info) header +=
" # [4] mean separation # [5] standard deviation of the separation distribution # [6] mean redshift # [7] standard deviation of the redshift distribution";
234 data_SS[i]->write(dir_output_resample, file, header, 10, 0);
237 vector<double> dd; data_SS[i]->get_data(dd);
239 xi_SubSamples.push_back(dd);
244 if (estimator==Estimator::_natural_)
245 m_dataset = correlation_NaturalEstimator(m_dd, m_rr);
246 else if (estimator==Estimator::_LandySzalay_)
247 m_dataset = correlation_LandySzalayEstimator(m_dd, m_rr, m_dr);
249 ErrorCBL(
"the chosen estimator is not implemented!",
"measureBootstrap",
"TwoPointCorrelation1D_angular.cpp");
251 m_dataset->set_covariance(covariance);
The class TwoPointCorrelation1D_angular.
void measureBootstrap(const int nMocks, const std::string dir_output_pairs=par::defaultString, const std::vector< std::string > dir_input_pairs={}, const std::string dir_output_resample=par::defaultString, const bool count_dd=true, const bool count_rr=true, const bool count_dr=true, const bool tcount=true, const Estimator estimator=Estimator::_LandySzalay_, const double fact=0.1, const int seed=3213) override
measure the angular two-point correlation function estimating the covariance with Bootstrap resamplin...
void measurePoisson(const std::string dir_output_pairs=par::defaultString, const std::vector< std::string > dir_input_pairs={}, const bool count_dd=true, const bool count_rr=true, const bool count_dr=true, const bool tcount=true, const Estimator estimator=Estimator::_LandySzalay_, const double fact=0.1) override
measure the angular two-point correlation function with Poisson errors
void measureJackknife(const std::string dir_output_pairs=par::defaultString, const std::vector< std::string > dir_input_pairs={}, const std::string dir_output_resample=par::defaultString, const bool count_dd=true, const bool count_rr=true, const bool count_dr=true, const bool tcount=true, const Estimator estimator=Estimator::_LandySzalay_, const double fact=0.1) override
measure the angular two-point correlation function estimating the covariance with Jackknife resamplin...
void read(const std::string dir, const std::string file) override
read the angular two-point correlation function
void measure(const ErrorType errorType=ErrorType::_Poisson_, const std::string dir_output_pairs=par::defaultString, const std::vector< std::string > dir_input_pairs={}, const std::string dir_output_resample=par::defaultString, const int nMocks=0, const bool count_dd=true, const bool count_rr=true, const bool count_dr=true, const bool tcount=true, const Estimator estimator=Estimator::_LandySzalay_, const double fact=0.1, const int seed=3213) override
measure the angular two-point correlation function
void set_parameters(const BinType binType, const double thetaMin, const double thetaMax, const int nbins, const double shift, const CoordinateUnits angularUnits=CoordinateUnits::_radians_, std::function< double(double)> angularWeight=nullptr, const bool compute_extra_info=false)
set the binning parameters
void write(const std::string dir=par::defaultString, const std::string file=par::defaultString, const int rank=0) const override
write the angular two-point correlation function
static const char fINT[]
conversion symbol for: int -> std::string
static const std::string defaultString
default std::string value
Estimator
the two-point correlation estimator
ErrorType
the two-point correlation function error type
The global namespace of the CosmoBolognaLib
std::string conv(const T val, const char *fact)
convert a number to a std::string
void checkDim(const std::vector< T > vect, const int val, const std::string vector, bool equal=true)
check if the dimension of a std::vector is equal/lower than an input value
int ErrorCBL(const std::string msg, const std::string functionCBL, const std::string fileCBL, const cbl::glob::ExitCode exitCode=cbl::glob::ExitCode::_error_)
throw an exception: it is used for handling exceptions inside the CosmoBolognaLib
void covariance_matrix(const std::vector< std::vector< double >> mat, std::vector< std::vector< double >> &cov, const bool JK=false)
compute the covariance matrix from an input dataset
CoordinateUnits
the coordinate units