DAQ/TPCRawDataTree/TPCRawDataTree.cc

0001
0002 #include "TPCRawDataTree.h"
0003
0004 #include <fun4all/Fun4AllReturnCodes.h>
0005 #include <phool/PHCompositeNode.h>
0006 #include <phool/PHIODataNode.h>    // for PHIODataNode
0007 #include <phool/PHNodeIterator.h>  // for PHNodeIterator
0008 #include <phool/PHObject.h>        // for PHObject
0009 #include <phool/getClass.h>
0010 #include <phool/phool.h>
0011
0012 #include <Event/Event.h>
0013 #include <Event/EventTypes.h>
0014 #include <Event/packet.h>
0015
0016 #include <TFile.h>
0017 #include <TTree.h>
0018 #include <TH1F.h>
0019 #include <TH2F.h>
0020
0021 #include <cassert>
0022 #include <iostream>
0023 #include <memory>
0024
0025 //____________________________________________________________________________..
0026 TPCRawDataTree::TPCRawDataTree(const std::string &name)
0027   : SubsysReco("TPCRawDataTree")
0028   , m_fname(name)
0029 {
0030   // reserve memory for max ADC samples
0031   m_adcSamples.resize(1024, 0);
0032   M.setMapNames("AutoPad-R1-RevA.sch.ChannelMapping.csv", "AutoPad-R2-RevA-Pads.sch.ChannelMapping.csv", "AutoPad-R3-RevA.sch.ChannelMapping.csv");
0033 }
0034
0035 //____________________________________________________________________________..
0036 int TPCRawDataTree::InitRun(PHCompositeNode *)
0037 {
0038   sectorNum = m_fname;
0039   size_t pos = sectorNum.find("TPC_ebdc");
0040   sectorNum.erase(sectorNum.begin(),sectorNum.begin()+pos+8);
0041   sectorNum.erase(sectorNum.begin()+2,sectorNum.end());
0042   if(sectorNum.at(0) == '0') sectorNum.erase(sectorNum.begin(),sectorNum.begin()+1);
0043   if(stoi(sectorNum) > 11) side = 1;
0044
0045   m_file = TFile::Open(m_fname.c_str(), "recreate");
0046   assert(m_file->IsOpen());
0047
0048   m_PacketTree = new TTree("PacketTree", "Each entry is one packet");
0049
0050   m_PacketTree->Branch("packet", &m_packet, "packet/I");
0051   m_PacketTree->Branch("frame", &m_frame, "frame/I");
0052   m_PacketTree->Branch("nWaveormInFrame", &m_nWaveormInFrame, "nWaveormInFrame/I");
0053   m_PacketTree->Branch("nTaggerInFrame", &m_nTaggerInFrame, "nTaggerInFrame/I");
0054   m_PacketTree->Branch("maxFEECount", &m_maxFEECount, "maxFEECount/I");
0055
0056   m_SampleTree = new TTree("SampleTree", "Each entry is one waveform");
0057
0058   m_SampleTree->Branch("packet", &m_packet, "packet/I");
0059   m_SampleTree->Branch("frame", &m_frame, "frame/I");
0060   m_SampleTree->Branch("nWaveormInFrame", &m_nWaveormInFrame, "nWaveormInFrame/I");
0061   m_SampleTree->Branch("maxFEECount", &m_maxFEECount, "maxFEECount/I");
0062   m_SampleTree->Branch("nSamples", &m_nSamples, "nSamples/I");
0063   m_SampleTree->Branch("adcSamples", &m_adcSamples[0], "adcSamples[nSamples]/s");
0064   m_SampleTree->Branch("fee", &m_fee, "fee/I");
0065   m_SampleTree->Branch("sampaAddress", &m_sampaAddress, "sampaAddress/I");
0066   m_SampleTree->Branch("sampaChannel", &m_sampaChannel, "sampaChannel/I");
0067   m_SampleTree->Branch("Channel", &m_Channel, "Channel/I");
0068   m_SampleTree->Branch("BCO", &m_BCO, "BCO/I");
0069   m_SampleTree->Branch("checksum", &m_checksum, "checksum/I");
0070   m_SampleTree->Branch("checksumError", &m_checksumError, "checksumError/I");
0071   m_SampleTree->Branch("Lv1TaggerBCO", &m_Lv1TaggerBCO, "Lv1TaggerBCO/l");
0072
0073   m_TaggerTree = new TTree("TaggerTree", "Each entry is one tagger for level 1 trigger or endat tag");
0074
0075   m_TaggerTree->Branch("packet", &m_packet, "packet/I");
0076   m_TaggerTree->Branch("frame", &m_frame, "frame/I");
0077   m_TaggerTree->Branch("tagger_type", &m_tagger_type, "tagger_type/s");
0078   m_TaggerTree->Branch("is_endat", &m_is_endat, "is_endat/b");
0079   m_TaggerTree->Branch("is_lvl1", &m_is_lvl1, "is_lvl1/b");
0080   m_TaggerTree->Branch("bco", &m_bco, "bco/l");
0081   m_TaggerTree->Branch("lvl1_count", &m_lvl1_count, "lvl1_count/i");
0082   m_TaggerTree->Branch("endat_count", &m_endat_count, "endat_count/i");
0083   m_TaggerTree->Branch("last_bco", &m_last_bco, "last_bco/l");
0084   m_TaggerTree->Branch("modebits", &m_modebits, "modebits/b");
0085
0086   R1_hist = new TH1F("R1_hist","R1_hist",1024,-0.5,1023.5);
0087   R2_hist = new TH1F("R2_hist","R2_hist",1024,-0.5,1023.5);
0088   R3_hist = new TH1F("R3_hist","R3_hist",1024,-0.5,1023.5);
0089
0090   R1_time = new TH2F("R1_time","R1_time",360,-0.5,359.5,1024,-0.5,1023.5);
0091   R2_time = new TH2F("R2_time","R2_time",360,-0.5,359.5,1024,-0.5,1023.5);
0092   R3_time = new TH2F("R3_time","R3_time",360,-0.5,359.5,1024,-0.5,1023.5);
0093
0094   TotalFEE = new TH1F("TotalFEE", "Total FEE", 26, -0.5, 25.5);
0095   TotalFEEsampa = new TH1F("TotalFEEsampa", "Total FEE + sampa", 26*8, -0.5, 25*8-.5);
0096   TotalFRAME = new TH1F("TotalFRAME", "Total FRAME", 21, -0.5, 20.5);
0097
0098   checksumError_fee = new TH1F("FEEWithError", "FEE with Error", 26, -0.5, 25.5);
0099   checksumError_feesampa = new TH1F("FEEsampaWithError", "FEE*8+sampa with Error", 26*8, -0.5, 25*8-.5);
0100   checksumError_frame = new TH1F("FRAMEWithError", "FRAME with Error", 21, -0.5, 20.5);
0101
0102   if (m_includeXYPos)
0103   {
0104     m_SampleTree->Branch("xPos", &m_xPos, "xPos/d");
0105     m_SampleTree->Branch("yPos", &m_yPos, "yPos/d");
0106   }
0107
0108   return Fun4AllReturnCodes::EVENT_OK;
0109 }
0110
0111 //____________________________________________________________________________..
0112 int TPCRawDataTree::process_event(PHCompositeNode *topNode)
0113 {
0114   Event *_event = findNode::getClass<Event>(topNode, "PRDF");
0115   if (_event == nullptr)
0116   {
0117     std::cout << "TPCRawDataTree::Process_Event - Event not found" << std::endl;
0118     return -1;
0119   }
0120   if (_event->getEvtType() >= 8)  /// special events
0121   {
0122     return Fun4AllReturnCodes::DISCARDEVENT;
0123   }
0124
0125   m_frame = _event->getEvtSequence();
0126
0127   for (int packet : m_packets)
0128   {
0129     if (Verbosity())
0130     {
0131       std::cout << __PRETTY_FUNCTION__ << " : decoding packet " << packet << std::endl;
0132     }
0133
0134     m_packet = packet;
0135
0136     std::unique_ptr<Packet> p(_event->getPacket(m_packet));
0137     if (!p)
0138     {
0139       if (Verbosity())
0140       {
0141         std::cout << __PRETTY_FUNCTION__ << " : missing packet " << packet << std::endl;
0142       }
0143
0144       continue;
0145     }
0146
0147     m_nWaveormInFrame = p->iValue(0, "NR_WF");
0148     m_nTaggerInFrame = p->lValue(0, "N_TAGGER");
0149     m_maxFEECount = p->iValue(0, "MAX_FEECOUNT");
0150
0151     for (int t = 0; t < m_nTaggerInFrame; t++)
0152     {
0153       /*uint16_t*/ m_tagger_type = (uint16_t) (p->lValue(t, "TAGGER_TYPE"));
0154       /*uint8_t*/ m_is_endat = (uint8_t) (p->lValue(t, "IS_ENDAT"));
0155       /*uint8_t*/ m_is_lvl1 = (uint8_t) (p->lValue(t, "IS_LEVEL1_TRIGGER"));
0156       /*uint64_t*/ m_bco = (uint64_t) (p->lValue(t, "BCO"));
0157       /*uint32_t*/ m_lvl1_count = (uint32_t) (p->lValue(t, "LEVEL1_COUNT"));
0158       /*uint32_t*/ m_endat_count = (uint32_t) (p->lValue(t, "ENDAT_COUNT"));
0159       /*uint64_t*/ m_last_bco = (uint64_t) (p->lValue(t, "LAST_BCO"));
0160       /*uint8_t*/ m_modebits = (uint8_t) (p->lValue(t, "MODEBITS"));
0161
0162       m_TaggerTree->Fill();
0163
0164       if (m_is_lvl1)
0165       {
0166         m_Lv1TaggerBCO = m_bco;
0167         m_Lv1TaggerCount = m_lvl1_count;
0168       }
0169     }
0170
0171     for (int wf = 0; wf < m_nWaveormInFrame; wf++)
0172     {
0173       m_BCO = p->iValue(wf, "BCO");
0174       m_nSamples = p->iValue(wf, "SAMPLES");
0175       m_fee = p->iValue(wf, "FEE");
0176
0177       m_sampaAddress = p->iValue(wf, "SAMPAADDRESS");
0178       m_sampaChannel = p->iValue(wf, "SAMPACHANNEL");
0179       m_Channel = p->iValue(wf, "CHANNEL");
0180       m_checksum = p->iValue(wf, "CHECKSUM");
0181       m_checksumError = p->iValue(wf, "CHECKSUMERROR");
0182
0183       TH1F *fillHist;
0184       TH2F *fillHist2D;
0185
0186       if(m_fee == 2 ||
0187          m_fee == 4 ||
0188          m_fee == 3 ||
0189          m_fee == 13 ||
0190          m_fee == 17 ||
0191          m_fee == 16){ fillHist=R1_hist; fillHist2D=R1_time;}
0192       else if(m_fee == 11 ||
0193          m_fee == 12 ||
0194          m_fee == 19 ||
0195          m_fee == 18 ||
0196          m_fee == 01 ||
0197          m_fee == 00 ||
0198          m_fee == 16 ||
0199          m_fee == 15){ fillHist=R2_hist; fillHist2D=R2_time;}
0200       else{ fillHist=R3_hist; fillHist2D=R3_time;}
0201
0202
0203       assert(m_nSamples <= (int) m_adcSamples.size());  // no need for movements in memory allocation
0204       for (int s = 0; s < m_nSamples; s++)
0205       {
0206         m_adcSamples[s] = p->iValue(wf, s);
0207         if(m_checksumError==0){
0208            fillHist->Fill(m_adcSamples[s]);
0209            fillHist2D->Fill(s,m_adcSamples[s]);
0210         }
0211         else {
0212           checksumError_fee->Fill(m_fee);
0213           checksumError_feesampa->Fill((m_fee*8. + m_sampaAddress));
0214           checksumError_frame->Fill(m_frame);
0215         }
0216         TotalFEE->Fill(m_fee);
0217         TotalFEEsampa->Fill((m_fee*8. + m_sampaAddress));
0218         TotalFRAME->Fill(m_frame);
0219       }
0220       if(m_includeXYPos)
0221       {
0222         int feeM = FEE_map[m_fee];
0223         if(FEE_R[m_fee]==2) feeM += 6;
0224         if(FEE_R[m_fee]==3) feeM += 14;
0225         int layer = M.getLayer(feeM, m_Channel);
0226         if(layer!=0)
0227         {
0228           double R = M.getR(feeM, m_Channel);
0229           double phi = (stoi(sectorNum)<12? -M.getPhi(feeM, m_Channel) : M.getPhi(feeM, m_Channel) ) + (stoi(sectorNum) - side*12. )* M_PI / 6. ;
0230           R /= 10.; //convert mm to cm
0231
0232           m_xPos = R*cos(phi);
0233           m_yPos = R*sin(phi);
0234         }
0235         else
0236         {
0237           m_xPos = 0.;
0238           m_yPos = 0.;
0239         }
0240       }
0241       m_SampleTree->Fill();
0242     }
0243
0244     m_PacketTree->Fill();
0245   }  //   for (int packet : m_packets)
0246
0247   return Fun4AllReturnCodes::EVENT_OK;
0248 }
0249
0250 //____________________________________________________________________________..
0251 int TPCRawDataTree::End(PHCompositeNode * /*topNode*/)
0252 {
0253   checksumError_fee->Divide(TotalFEE);
0254   checksumError_feesampa->Divide(TotalFEEsampa);
0255   checksumError_frame->Divide(TotalFRAME);
0256
0257   TotalFEE->SetDirectory(0);
0258   TotalFEEsampa->SetDirectory(0);
0259   TotalFRAME->SetDirectory(0);
0260
0261   m_file->Write();
0262
0263   std::cout << __PRETTY_FUNCTION__ << " : completed saving to " << m_file->GetName() << std::endl;
0264   m_file->ls();
0265
0266   m_file->Close();
0267
0268   return Fun4AllReturnCodes::EVENT_OK;
0269 }