calibrations/mbd/read_dstmbd.C

0001 #ifndef MACRO_READ_DSTMBD_C
0002 #define MACRO_READ_DSTMBD_C
0003 //
0004 // macro to read in MBD data after the waveforms are processed
0005 // at this stage, we only have charge and time from each channel, not the full waveform
0006 //
0007 #include <mbd/MbdDefs.h>
0008 #include <mbd/MbdOut.h>
0009 #include <mbd/MbdPmtContainer.h>
0010 #include <mbd/MbdPmtHit.h>
0011
0012 #include <TChain.h>
0013 #include <TFile.h>
0014
0015 #include <fstream>
0016
0017 R__LOAD_LIBRARY(libmbd_io.so)
0018
0019 // Set up variables to read from TTree
0020 Int_t   f_evt{0};
0021 UShort_t f_clk{0};
0022 UShort_t f_femclk{0};
0023 Short_t  f_npmt;
0024 Float_t f_tt[MbdDefs::MBD_N_PMT];  // time from t-channels
0025 Float_t f_tq[MbdDefs::MBD_N_PMT];  // time from q-channels
0026 Float_t f_q[MbdDefs::MBD_N_PMT];   // charge
0027
0028 Short_t f_bn[MbdDefs::MBD_N_ARMS]; // num hit PMTs
0029 Float_t f_bq[MbdDefs::MBD_N_ARMS]; // chargesum
0030 Float_t f_bt[MbdDefs::MBD_N_ARMS]; // mean time in arm
0031 Float_t f_bz;                      // z-vertex
0032 Float_t f_bt0;                     // t-zero
0033
0034 MbdOut *mbdout{nullptr};
0035 MbdPmtContainer *mbdpmts{nullptr};
0036
0037 TFile *tfile {nullptr};
0038 TTree *tree {nullptr};
0039
0040 void dstmbd_GetEntry(const int ientry)
0041 {
0042     tree->GetEntry(ientry);
0043
0044     f_evt = mbdout->get_evt();
0045     f_clk = mbdout->get_clock();
0046     f_femclk = mbdout->get_femclock();
0047     f_bz = mbdout->get_zvtx();
0048     f_bt0 = mbdout->get_t0();
0049     for (int iarm=0; iarm<2; iarm++)
0050     {
0051       f_bn[iarm] = mbdout->get_npmt(iarm);
0052       f_bq[iarm] = mbdout->get_q(iarm);
0053       f_bt[iarm] = mbdout->get_time(iarm);
0054     }
0055     f_npmt = f_bn[0] + f_bn[1];
0056
0057     // clear pmt info
0058     for (int ipmt=0; ipmt<128; ipmt++)
0059     {
0060       f_q[ipmt] = 0.;
0061       f_tt[ipmt] = NAN;
0062       f_tq[ipmt] = NAN;
0063     }
0064
0065     // Loop over each MBD/BBC PMT container
0066     for (int ipmt=0; ipmt<mbdpmts->get_npmt(); ipmt++)
0067     {
0068       MbdPmtHit* mbdpmt = mbdpmts->get_pmt(ipmt);
0069
0070       f_q[ipmt] = mbdpmt->get_q();
0071       f_tt[ipmt] = mbdpmt->get_tt();
0072       f_tq[ipmt] = mbdpmt->get_tq();
0073     }
0074 }
0075
0076 void read_dstmbd(const char *tfname = "beam_mbd-00009184-0000_mbd.root")
0077 {
0078   std::cout << "tfname " << tfname << std::endl;
0079
0080   // Set up TTree
0081 //  int is_dst = 0; // whether reading from DST or private root files
0082   tfile = new TFile(tfname,"READ");
0083   tree = (TTree*)tfile->Get("T");
0084
0085   tree->SetBranchAddress("DST#MBD#MbdOut",&mbdout);
0086   tree->SetBranchAddress("DST#MBD#MbdPmtContainer",&mbdpmts);
0087
0088   /* older dst style, with better leaf level
0089   tree->SetBranchAddress("EvtSequence",&f_evt);
0090   tree->SetBranchAddress("clk", &f_clk);
0091   tree->SetBranchAddress("femclk", &f_femclk);
0092   tree->SetBranchAddress("bns",&f_bn[0]);
0093   tree->SetBranchAddress("bnn",&f_bn[1]);
0094   tree->SetBranchAddress("bqs",&f_bq[0]);
0095   tree->SetBranchAddress("bqn",&f_bq[1]);
0096   tree->SetBranchAddress("bts",&f_bt[0]);
0097   tree->SetBranchAddress("btn",&f_bt[1]);
0098   tree->SetBranchAddress("bz",&f_bz);
0099   tree->SetBranchAddress("bt0",&f_bt0);
0100   tree->SetBranchAddress("MbdPmtHits.bq",f_q);
0101   tree->SetBranchAddress("MbdPmtHits.btt",f_tt);
0102   tree->SetBranchAddress("MbdPmtHits.btq",f_tq);
0103   */
0104
0105   // Event loop, each ientry is one triggered event
0106   int nentries = tree->GetEntries();
0107   std::cout << "checking output of dst file " << tfname << std::endl;
0108   std::cout << "nentries = " << nentries << std::endl;
0109   for (int ientry=0; ientry<10; ientry++)
0110   {
0111     //std::cout << "ientry " << ientry << std::endl;
0112     dstmbd_GetEntry(ientry);
0113     //tree->GetEntry(ientry+100);
0114
0115     if (ientry<10)
0116     {
0117       // print charge from channels 0 and 127
0118       std::cout << f_evt << "\t" << f_q[0] << "\t" << f_tt[0] << std::endl;
0119       std::cout <<  "\t" << f_q[127] << "\t" << f_tt[127] << std::endl;
0120       std::cout <<  "\t" << f_bz << std::endl;
0121     }
0122
0123   }
0124
0125 }
0126 #endif