Skip to content
Snippets Groups Projects
Commit 8b410851 authored by Chao Peng's avatar Chao Peng
Browse files

Follow up: remove an old script and fix some minor issues

parent 71d428f8
Branches
No related tags found
No related merge requests found
......@@ -4,7 +4,7 @@ import getpass
import paramiko
root_path = os.path.join(os.path.dirname(os.path.realpath(__file__)), '..')
root_path = os.path.join(os.path.dirname(os.path.realpath(__file__)), '../..')
exec_path = os.path.join(root_path, 'bin/analyze')
......@@ -21,10 +21,10 @@ timing = os.path.join(root_path, 'config/mapmt_timing.conf')
run_options = [
'--config-module=%s' % modules,
'--config-timing=%s' % timing,
# '--config-timing=%s' % timing,
# '-c', str(tcut),
# '-f', str(fire_thres),
'--raw-output'
# '--raw-output'
]
# runs
......
......@@ -73,7 +73,7 @@ int main(int argc, char* argv[])
}
int nev = -1;
std::string mconf = "config/esb_module.conf";
std::string mconf = "config/mapmt_module.conf";
for (auto &opt : conf_opt.GetOptions()) {
switch (opt.mark) {
......@@ -231,7 +231,7 @@ void fill_branch(const Fadc250Data &slot_data, const Module &mod, std::unordered
// fill branch data into the root tree
void fill_tree(TTree *tree, std::unordered_map<std::string, BranchData> &brdata, bool &init, int mode)
{
if ((mode != 1) && (mode != 3)) {
if ((mode != 1) && (mode != 3) && (mode != 10)) {
std::cout << "Warning: unsupported mode " << mode << ", data won't be recorded." << std::endl;
return;
}
......@@ -299,8 +299,8 @@ void write_raw_data(const std::string &dpath, const std::string &opath, int nev,
if ((evdata->GetNumRaw(mod.crate, mod.slot) == 0) || (fadc == nullptr)) { continue; }
auto slot_data = pack_fadc250_data(fadc);
// check data mode
if (slot_data.GetMode() > 0 && data_mode < 0) {
data_mode = slot_data.GetMode();
if (slot_data.mode > 0 && data_mode < 0) {
data_mode = slot_data.mode;
}
// fill branch data
fill_branch(slot_data, mod, brdata);
......
......@@ -125,7 +125,7 @@ int main(int argc, char* argv[])
conf_opt.SetDesc("usage: %0 <data_file> <out_file>");
conf_opt.SetDesc('n', "number of events to process (< 0 means all).");
conf_opt.SetDesc('m', "configuration file for modules to be read-in, default \"config/mapmt_module.conf\"");
conf_opt.SetDesc('m', "configuration file for modules to be read-in, default \"config/esb_module.conf\"");
if (!conf_opt.ParseArgs(argc, argv) || conf_opt.NbofArgs() != 2) {
std::cout << conf_opt.GetInstruction() << std::endl;
......@@ -342,6 +342,7 @@ void write_raw_data(const std::string &dpath, const std::string &opath, int nev,
int count = 0;
bool init_tree = false;
int data_mode = -1;
while ((datafile.codaRead() == S_SUCCESS) && (nev-- != 0)) {
// read-in data
evdata->LoadEvent(datafile.getEvBuffer());
......@@ -359,7 +360,6 @@ void write_raw_data(const std::string &dpath, const std::string &opath, int nev,
int blvl = 1;
simpleGetRocBlockLevel(modules.front().crate, FADC_BANK, &blvl);
int data_mode = -1;
for (int ii = 0; ii < blvl; ++ii) {
// clear data buffer
for (auto &br : brdata) {
......@@ -387,7 +387,7 @@ void write_raw_data(const std::string &dpath, const std::string &opath, int nev,
auto slot_data = fadc250_decoder(header, dbuf, len);
// check data mode
if (slot_data.GetMode() > 0 && data_mode < 0) {
if (data_mode < 0 && slot_data.GetMode() > 0) {
data_mode = slot_data.GetMode();
}
// fill branch data
......
#include "utils.h"
#include "ConfigParser.h"
#include "ConfigObject.h"
#include "ConfigOption.h"
#include "THaCodaFile.h"
#include "THaEvData.h"
#include "CodaDecoder.h"
#include "Fadc250Module.h"
#include "evio.h"
#include "TTree.h"
#include "TFile.h"
#include "TH1.h"
#define PROGRESS_COUNT 1000
std::vector<Module> read_modules(const std::string &path);
std::unordered_map<std::string, std::vector<float>> read_timing(const std::string &path);
void process_data(const std::string &dpath, const std::string &opath, int nev, const std::vector<Module> &modules,
std::unordered_map<std::string, std::vector<float>> timing, float nsig, int nhit);
int main(int argc, char* argv[])
{
// setup input arguments
ConfigOption conf_opt;
conf_opt.AddOpts(ConfigOption::help_message, 'h', "help");
conf_opt.AddOpt(ConfigOption::arg_require, 'n');
conf_opt.AddOpt(ConfigOption::arg_require, 'c');
conf_opt.AddOpt(ConfigOption::arg_require, 'i');
conf_opt.AddLongOpt(ConfigOption::arg_require, "config-module", 'm');
conf_opt.AddLongOpt(ConfigOption::arg_require, "config-timing", 't');
conf_opt.SetDesc("usage: %0 <data_file> <out_dir>");
conf_opt.SetDesc('n', "number of events to process (< 0 means all).");
conf_opt.SetDesc('i', "threshold on number of fired modules, default 3.");
conf_opt.SetDesc('c', "timing cut factor, default 3 (3 sigmas).");
conf_opt.SetDesc('m', "configuration file for modules to be read-in, default \"config/mapmt_module.conf\"");
conf_opt.SetDesc('t', "configuration file for signal timing to be read-in, default \"config/mapmt_timing.conf\"");
if (!conf_opt.ParseArgs(argc, argv) || conf_opt.NbofArgs() != 2) {
std::cout << conf_opt.GetInstruction() << std::endl;
return -1;
}
int nev = -1;
int nhit = 3;
float nsig = 3.;
std::string mconf = "config/mapmt_module.conf";
std::string tconf = "config/mapmt_timing.conf";
for (auto &opt : conf_opt.GetOptions()) {
switch (opt.mark) {
case 'n':
nev = opt.var.Int();
break;
case 'c':
nsig = opt.var.Float();
break;
case 'i':
nhit = opt.var.Int();
break;
case 'm':
mconf = opt.var.String();
break;
case 't':
tconf = opt.var.String();
break;
default :
std::cout << conf_opt.GetInstruction() << std::endl;
return -1;
}
}
auto modules = read_modules(mconf);
auto timing = read_timing(tconf);
process_data(conf_opt.GetArgument(0), conf_opt.GetArgument(1), nev, modules, timing, nsig, nhit);
}
std::vector<Module> read_modules(const std::string &path)
{
// read in file
std::string buffer = ConfigParser::file_to_string(path);
// remove comments
ConfigParser::comment_line(buffer, "#", "\n");
// get content blocks
auto text = ConfigParser::break_into_blocks(buffer, "{", "}");
ConfigObject conf;
std::vector<Module> res;
for(auto &b : text.blocks)
{
// module
if (!ConfigParser::case_ins_equal("Module", b.label))
continue;
auto btext = ConfigParser::break_into_blocks(b.content, "{", "}");
conf.ReadConfigString(btext.residual);
// module attributes
Module mod;
mod.crate = conf.GetConfigValue("crate").Int();
mod.slot = conf.GetConfigValue("slot").Int();
mod.type = str2ModuleType(conf.GetConfigValue("type").c_str());
// channels
for (auto &sub : btext.blocks) {
if (!ConfigParser::case_ins_equal("Channels", sub.label))
continue;
ConfigParser parser;
parser.ReadBuffer(sub.content.c_str());
while(parser.ParseLine()) {
mod.channels.emplace_back(Channel{
parser.TakeFirst().Int(),
parser.TakeFirst().String(),
str2ChannelType(parser.TakeFirst().c_str())
});
}
}
res.emplace_back(mod);
}
// print out all channels
std::cout << "Read-in " << res.size() << " modules from \"" << path << "\"." << std::endl;
for (auto &mod : res) {
std::cout << "Module: crate = " << mod.crate
<< ", slot = " << mod.slot
<< ", type = " << ModuleType2str(mod.type)
<< std::endl;
for (auto &ch : mod.channels) {
std::cout << "\t Channel " << ch.id
<< ": name = " << ch.name
<< ", type = " << ChannelType2str(ch.type)
<< std::endl;
}
}
return res;
}
std::unordered_map<std::string, std::vector<float>> read_timing(const std::string &path)
{
std::unordered_map<std::string, std::vector<float>> res;
ConfigParser parser;
std::cout << "Reading channel timing info from \"" << path << "\"." << std::endl;
parser.ReadFile(path);
std::string name;
float center, sigma;
while (parser.ParseLine()) {
parser >> name >> center >> sigma;
if (res.find(name) != res.end()) {
std::cout << "WARNING: Duplicated timing information for channel " << name << std::endl;
}
res[name] = {center, sigma};
std::cout << name << ": (" << center << ", " << sigma << ")" << std::endl;
}
return res;
}
bool has_sample(const std::unordered_map<std::string, std::vector<uint32_t>> &samples)
{
for (auto &sample : samples) {
if (sample.second.size() > 0) return true;
}
return false;
}
void process_data(const std::string &dpath, const std::string &opath, int nev, const std::vector<Module> &modules,
std::unordered_map<std::string, std::vector<float>> timing, float nsig, int nhit)
{
Decoder::THaCodaFile datafile(dpath.c_str());
THaEvData *evdata = new Decoder::CodaDecoder();
std::unordered_map<std::string, std::vector<PulseData>> event;
std::unordered_map<std::string, std::vector<uint32_t>> samples;
std::vector<std::string> signals;
for (auto &mod : modules) {
for (auto &ch : mod.channels) {
if (ch.type == kSignal)
signals.push_back(ch.name);
}
}
std::ofstream output(opath);
char csv_sep = ',';
// csv file header
output << "Event Number" << csv_sep << "Fired Counts" << csv_sep <<"Signal Sum" << csv_sep << "LED";
for (auto &sig : signals) {
output << csv_sep << sig;
}
output << std::endl;
int count = 0;
float int_thres = 0.;
while ((datafile.codaRead() == S_SUCCESS) && (nev-- != 0)) {
// read-in data
evdata->LoadEvent(datafile.getEvBuffer());
if((++count % PROGRESS_COUNT) == 0) {
std::cout << "Processed events - " << count << "\r" << std::flush;
}
// read data into event buffer
for (auto &mod : modules) {
auto *fadc = dynamic_cast<Decoder::Fadc250Module*>(evdata->GetModule(mod.crate, mod.slot));
bool empty_slot = ((evdata->GetNumRaw(mod.crate, mod.slot) == 0) || (fadc == nullptr));
for (auto &ch : mod.channels) {
if (empty_slot) {
event[ch.name] = std::vector<PulseData>();
samples[ch.name] = std::vector<uint32_t>();
} else {
event[ch.name] = get_pulses(fadc, ch.id, mod.type);
if (fadc->GetNumSamples(ch.id) > 0) {
samples[ch.name] = fadc->GetPulseSamplesVector(ch.id);
} else {
samples[ch.name] = std::vector<uint32_t>();
}
}
}
}
// further process event
// no reference time
if (event["Ref"].size() == 0)
continue;
auto ref_time = event["Ref"].front().time;
// Calosum timing
float trg_pulse = 0.;
float trg_time = 0.;
for (auto &cal : event["CaloSum"]) {
if (trg_pulse > cal.integral)
continue;
trg_time = cal.time;
trg_pulse = cal.integral;
}
if (trg_time == 0.)
continue;
int fire = 0;
float sum = 0;
for (auto &sig : signals) {
float integral = int_thres;
for (auto &ev : event[sig]) {
float tdiff = ev.time - trg_time;
if ((std::abs(tdiff - timing[sig][0]) < nsig*timing[sig][1]) &&
(ev.integral > integral)) {
integral = ev.integral;
}
}
if (integral > int_thres) {
fire ++;
sum += integral;
}
}
bool led_trigger = (event["LED"].size() > 0);
if (led_trigger || ((fire >= nhit) && has_sample(samples))) {
output << count << csv_sep << fire << csv_sep << sum << csv_sep << led_trigger;
for (auto &sig : signals) {
output << csv_sep << "\"[";
for (auto &spl : samples[sig]) {
output << spl << ",";
}
output << "]\"";
}
output << std::endl;
}
}
std::cout << "Processed events - " << count << std::endl;
output.close();
}
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment