#include #include #include #include #include #include #include #include "filval/filval.hpp" #include "filval/root/filval.hpp" #include #include "TrackingNtuple.h" #include "analysis/obj_types.cpp" using namespace std; using namespace fv; using namespace fv::root; using namespace std; #define HIT_TYPE_PIXEL 0 #define HIT_TYPE_GLUED 1 #define HIT_TYPE_STRIP 2 Value>* seeds; Value>* pixrec_hits; Value>* tracks; Value>* sim_hits; Value>* sim_tracks; void register_objects(TrackingDataSet& tds){ seeds = register_seeds(tds); pixrec_hits = register_pixrec_hits(tds); tracks = register_tracks(tds); sim_hits = register_sim_hits(tds); sim_tracks = register_sim_tracks(tds); } void setup_dphi_vs_eta(TrackingDataSet& tds){ auto& calc_dphi_vs_eta = func,vector>(vector, vector, vector, int)>("calc_dphi_vs_eta", FUNC(([](const vector& seeds, const vector pixrec_hits, const vector sim_hits, int layer){ vector dphi; vector eta; for(const Seed &seed : seeds){ int hitIdx, hitType; for(auto tup : boost::combine(seed.hitIdx, seed.hitType)){ boost::tie(hitIdx, hitType) = tup; if(hitType == HIT_TYPE_PIXEL){ const PixRecHit &hit = pixrec_hits[hitIdx]; if(layer && hit.lay != layer) continue; float hit_phi = atan2(hit.y, hit.x); for(const int& simHitIdx : hit.simHitIdx){ const SimHit &sim_hit = sim_hits[simHitIdx]; dphi.push_back(hit_phi-atan2(sim_hit.y, sim_hit.x)); eta.push_back(seed.eta); } } } } return std::make_pair(dphi,eta); }))); auto dphi_vs_eta_B1 = fv::apply(calc_dphi_vs_eta, fv::tuple(seeds, pixrec_hits, sim_hits, constant("1", 1)), "dphi_vs_eta_B1"); auto dphi_vs_eta_B2 = fv::apply(calc_dphi_vs_eta, fv::tuple(seeds, pixrec_hits, sim_hits, constant("2", 2)), "dphi_vs_eta_B2"); auto dphi_vs_eta_B3 = fv::apply(calc_dphi_vs_eta, fv::tuple(seeds, pixrec_hits, sim_hits, constant("3", 3)), "dphi_vs_eta_B3"); auto dphi_vs_eta_B4 = fv::apply(calc_dphi_vs_eta, fv::tuple(seeds, pixrec_hits, sim_hits, constant("4", 4)), "dphi_vs_eta_B4"); tds.register_container>("dphi_vs_eta_B1", "dphi vs eta (BPIX L1)", dphi_vs_eta_B1, 100, -0.002, 0.002, 20, -3, 3, "dphi", "eta"); tds.register_container>("dphi_vs_eta_B2", "dphi vs eta (BPIX L2)", dphi_vs_eta_B2, 100, -0.002, 0.002, 20, -3, 3, "dphi", "eta"); tds.register_container>("dphi_vs_eta_B3", "dphi vs eta (BPIX L3)", dphi_vs_eta_B3, 100, -0.002, 0.002, 20, -3, 3, "dphi", "eta"); tds.register_container>("dphi_vs_eta_B4", "dphi vs eta (BPIX L4)", dphi_vs_eta_B4, 100, -0.002, 0.002, 20, -3, 3, "dphi", "eta"); } void setup_dz_vs_eta(TrackingDataSet& tds){ auto& calc_dz_vs_eta = func,vector>(vector, vector, vector, int)>("calc_dz_vs_eta", FUNC(([](const vector& seeds, const vector pixrec_hits, const vector sim_hits, int layer){ vector dz; vector eta; for(const Seed &seed : seeds){ int hitIdx, hitType; for(auto tup : boost::combine(seed.hitIdx, seed.hitType)){ boost::tie(hitIdx, hitType) = tup; if(hitType == HIT_TYPE_PIXEL){ const PixRecHit &hit = pixrec_hits[hitIdx]; if(layer && hit.lay != layer) continue; for(const int& simHitIdx : hit.simHitIdx){ dz.push_back(hit.z-sim_hits[simHitIdx].z); eta.push_back(seed.eta); } } } } return std::make_pair(dz,eta); }))); auto dz_vs_eta_B1 = fv::apply(calc_dz_vs_eta, fv::tuple(seeds, pixrec_hits, sim_hits, constant("1", 1)), "dz_vs_eta_B1"); auto dz_vs_eta_B2 = fv::apply(calc_dz_vs_eta, fv::tuple(seeds, pixrec_hits, sim_hits, constant("2", 2)), "dz_vs_eta_B2"); auto dz_vs_eta_B3 = fv::apply(calc_dz_vs_eta, fv::tuple(seeds, pixrec_hits, sim_hits, constant("3", 3)), "dz_vs_eta_B3"); auto dz_vs_eta_B4 = fv::apply(calc_dz_vs_eta, fv::tuple(seeds, pixrec_hits, sim_hits, constant("4", 4)), "dz_vs_eta_B4"); tds.register_container>("dz_vs_eta_B1", "dz vs eta (BPIX L1)", dz_vs_eta_B1, 20, -0.01, .01, 20, -3, 3, "dz", "eta"); tds.register_container>("dz_vs_eta_B2", "dz vs eta (BPIX L2)", dz_vs_eta_B2, 20, -0.01, .01, 20, -3, 3, "dz", "eta"); tds.register_container>("dz_vs_eta_B3", "dz vs eta (BPIX L3)", dz_vs_eta_B3, 20, -0.01, .01, 20, -3, 3, "dz", "eta"); tds.register_container>("dz_vs_eta_B4", "dz vs eta (BPIX L4)", dz_vs_eta_B4, 20, -0.01, .01, 20, -3, 3, "dz", "eta"); } void run_analysis(const std::string& input_filename, const std::string& data_label, bool silent){ gSystem->Load("libfilval.so"); auto replace_suffix = [](const std::string& input, const std::string& new_suffix){ return input.substr(0, input.find_last_of(".")) + new_suffix; }; string log_filename = replace_suffix(input_filename, "_result.log"); fv::util::Log::init_logger(log_filename, fv::util::LogPriority::kLogDebug); string output_filename = replace_suffix(input_filename, "_result.root"); TrackingDataSet tds(output_filename, input_filename, data_label, "trackingNtuple/tree"); register_objects(tds); setup_dz_vs_eta(tds); setup_dphi_vs_eta(tds); tds.process(silent); tds.save_all(); } int main(int argc, char * argv[]) { fv::util::ArgParser args(argc, argv); if(!args.cmdOptionExists("-l") || !args.cmdOptionExists("-f")) { cout << "Usage: ./main (-s) -l DATA_LABEL -f treefile.root" << endl; return -1; } bool silent = args.cmdOptionExists("-s"); string input_filename = args.getCmdOption("-f"); string data_label = args.getCmdOption("-l"); run_analysis(input_filename, data_label, silent); }