#include "analyzer.h" #include "run.h" #include "bin.h" #include "histo.h" #include "const.h" #include #include #include #include #include using std::string; using namespace std; Analyzer::Analyzer() { reaction = UNDEF; list_of_runs = NULL; list_of_histograms = NULL; } Analyzer::~Analyzer() { } bool Analyzer::get_new_line (fstream& f, string& line) { // // reads the next non-empty and not-commeted-out line from file f // while (!f.eof()) { getline (f, line); unsigned int i = 0; while (i < line.length() && line.at(i) == ' ') i++; if (i < line.length() && line.at(i) != ' ' && line.at(i) != '#') return true; }; return false; }; // ---------------------------------------------------------------------------- bool Analyzer::read_input_file(string filename) { string line; fstream file(filename.c_str(), ios::in); if (filename.rfind(".palm") != string::npos) file_prefix = filename.substr(0, filename.rfind(".palm")); else file_prefix = filename; if (!file.is_open()) { cout << "Error reading palmetto input file: " << filename << endl; return false; } // Read reaction / kinematics flag if (get_new_line(file, line)) { // convert the line to upper case: std::transform(line.begin(), line.end(), line.begin(), (int(*)(int)) toupper); if (line.compare("EP_COINC") == 0) reaction = EP_COINC; else if (line.compare("EP_SINGLE") == 0) reaction = EP_SINGLE; else if (line.compare("EEP") == 0) reaction = EEP; else if (line.compare("GAMMAD") == 0 ) reaction = GAMMAD; else { cout << "ERROR: unknown reaction: " << line << endl; exit(1); } } // Read list of root files Run * last_run = NULL; bool skip = false; while (get_new_line(file, line) && line.compare("end") != 0) { if (line.compare("skip") == 0) skip = true; if (!skip) { Run *run_pointer = new Run(line); run_pointer->set_next_run(NULL); if (list_of_runs == NULL) list_of_runs = run_pointer; if (last_run != NULL) last_run->set_next_run(run_pointer); last_run = run_pointer; } } // cout << "List of runs ... " << endl; // cout << list_of_runs << endl; // Read beam parameters while (get_new_line(file, line) && line.compare("end") != 0) { int rmin, rmax; double h; sscanf (line.c_str(), "%d %d %lf", &rmin, &rmax, &h); Run * run_pointer = list_of_runs; while (run_pointer != NULL) { if (run_pointer->RunNumber >= rmin && run_pointer->RunNumber <= rmax) { run_pointer->beam_parameter_ok = true; run_pointer->h = h; } run_pointer = run_pointer->get_next_run(); } } // Read kinematics parameters while (get_new_line(file, line) && line.compare("end") != 0) { char name[20]; int rmin, rmax, skip; double e0, pLeft, thLeft, pRight, thRight, xfpp; char cstr[80]; string precession; sscanf (line.c_str(), "%s %d %d %lf %lf %lf %lf %lf %lf %s %d", name, &rmin, &rmax, &e0, &pLeft, &thLeft, &pRight, &thRight, &xfpp, cstr, &skip); precession = cstr; // convert precession to upper case: std::transform(precession.begin(), precession.end(), precession.begin(), (int(*)(int)) toupper); Run * run_pointer = list_of_runs; while (run_pointer != NULL) { if (run_pointer->RunNumber >= rmin && run_pointer->RunNumber <= rmax) { run_pointer->run_parameter_ok = true; run_pointer->kinematics = name; run_pointer->e0 = e0; run_pointer->pLeft = pLeft; // run_pointer->thLeft = thLeft * TMath::Pi() / 180.; run_pointer->thLeft = thLeft * 3.14159265 / 180.; run_pointer->pRight = pRight; // run_pointer->thRight = thRight * TMath::Pi() / 180.; run_pointer->thRight = thRight * 3.14159265 / 180.; run_pointer->xfpp = xfpp; run_pointer->skip = (skip == 1); if (precession == "DIPOLE") run_pointer->spin_transport_model = DIPOLE; else if (precession == "PENTCHEV") run_pointer->spin_transport_model = PENTCHEV; else if (precession == "PENTCHEV1") run_pointer->spin_transport_model = PENTCHEV1; else if (precession == "PENTCHEV2") run_pointer->spin_transport_model = PENTCHEV2; else if (precession.substr(precession.length()-5, precession.length()) == ".SPIN") { run_pointer->spin_transport_model = COSY; run_pointer->spin_transport_file = cstr; } else { cout << "ERROR: unknown spin transport model: " << precession << endl; exit(-1); } } run_pointer = run_pointer->get_next_run(); } } // Read cuts global_cut = ""; while (get_new_line(file, line) && line.compare("end") != 0) { TCut cut = line.c_str(); global_cut += cut; } // read histograms Histo * last_histo = NULL; while (get_new_line(file, line) && line.compare("end") != 0) { Histo *histo_pointer = new Histo(line); histo_pointer->set_next_histo(NULL); if (list_of_histograms == NULL) list_of_histograms = histo_pointer; if (last_histo != NULL) last_histo->set_next_histo(histo_pointer); last_histo = histo_pointer; } // check that all runs have proper information Run * run_pointer = list_of_runs; while (run_pointer != NULL) { if (!run_pointer->run_parameter_ok) { cout << "Error: missing run paramter for run: " << run_pointer->rootfile << endl; exit(1); } else if (!run_pointer->beam_parameter_ok) { cout << "Error: missing beam paramter for run: " << run_pointer->rootfile << endl; exit(1); } run_pointer = run_pointer->get_next_run(); } file.close(); return true; } // ---------------------------------------------------------------------------- void Analyzer::run_loop(string mode) { Run * run = list_of_runs; bool first_run_done = false; while (run != NULL) { if (!run->skip) { run->open(global_cut); list_of_histograms->init_all(run); while ( run->is_next_event() ) { run->calc_kinematics(reaction); run->calc_spin_transport(reaction); run->calc_weights(); list_of_histograms->add_all(run); summary.add(0, run); } run->close(); first_run_done = true; } run = run->get_next_run(); if (mode.compare("test") == 0 && first_run_done) break; } } // ---------------------------------------------------------------------------- const char* Analyzer::get_filename(string suffix) { return (file_prefix + suffix).c_str(); } // ---------------------------------------------------------------------------- void Analyzer::histogram_output() { Histo* h = list_of_histograms; while (h != NULL) { ofstream output; output.open (get_filename("." + h->get_type() + "." + h->get_name() + ".dat"), ios::out); output << *h; output.close(); h = h->get_next_histo(); } } // ---------------------------------------------------------------------------- ostream& operator<< (ostream &os, const Analyzer &obj) { time_t rawtime; struct tm * timeinfo; time (&rawtime); timeinfo = localtime (&rawtime); os << endl << "=====================================" << endl << endl << "P A L M E T T O - S u m m a r y" << endl << endl << "=====================================" << endl << endl; os << "Date: " << asctime (timeinfo) << endl; os << "Reaction: "; switch (obj.reaction) { case EP_SINGLE: os << "EP_SINGLE"; break; case EP_COINC: os << "EP_COINC"; break; case EEP: os << "EEP"; break; case GAMMAD : os << "GAMMAD"; break; default: os << "UNKNOWN"; break; }; os << endl << endl; const Run* a = obj.list_of_runs; os << a; os << endl << "Results:" << endl << "-------" << endl; Bin s = obj.summary; os << "Number of events = " << s.get_n() << endl; os << " Px = " << s.get_Px() << " +/- " << s.get_dPx() << endl; os << " Py = " << s.get_Py() << " +/- " << s.get_dPy() << endl; os << " Pz = " << s.get_Pz() << " +/- " << s.get_dPz() << endl; os << " R = " << s.get_R() << " +/- " << s.get_dR() << endl; os << "Analyzing power:" << endl; os << " model: Ac = " << s.get_Ac_calc() << endl; os << " data: Ac = " << s.get_Ac() << " +/- " << s.get_dAc() << endl; os << endl << "Histograms:" << endl << "----------" << endl; Histo* h = obj.list_of_histograms; while (h != NULL) { os << *h; h = h->get_next_histo(); }; return os; }