Belle II Software  release-05-02-19
CDCDQMModule.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2010 - Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Makoto Uchida *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 
11 // Own include
12 #include <cdc/modules/cdcDQM/CDCDQMModule.h>
13 
14 // CDC
15 
16 // framework - DataStore
17 #include <framework/datastore/StoreArray.h>
18 
19 // Dataobject classes
20 #include <framework/database/DBObjPtr.h>
21 
22 #include <TF1.h>
23 #include <TVector3.h>
24 #include <TDirectory.h>
25 
26 #include <fstream>
27 #include <math.h>
28 
29 using namespace std;
30 
31 namespace Belle2 {
37  //-----------------------------------------------------------------
38  // Register module
39  //-----------------------------------------------------------------
40 
41  REG_MODULE(CDCDQM);
42 
43  CDCDQMModule::CDCDQMModule() : HistoModule()
44  {
45  // set module description (e.g. insert text)
46  setDescription("Make summary of data quality.");
47  addParam("MinHits", m_minHits, "Include only events with more than MinHits hits in ARICH", 0);
49  }
50 
52  {
53  }
54 
56  {
57 
58  TDirectory* oldDir = gDirectory;
59 
60  oldDir->mkdir("CDC");
61  oldDir->cd("CDC");
62  m_hNEvents = new TH1F("hNEvents", "hNEvents", 10, 0, 10);
63  m_hNEvents->GetXaxis()->SetBinLabel(1, "number of events");
64  m_hOcc = new TH1F("hOcc", "hOccupancy", 150, 0, 1.5);
65  m_hADC = new TH2F("hADC", "hADC", 300, 0, 300, 1000, 0, 1000);
66  m_hADCTOTCut = new TH2F("hADCTOTCut", "hADCTOTCut", 300, 0, 300, 1000, 0, 1000);
67  m_hTDC = new TH2F("hTDC", "hTDC", 300, 0, 300, 1000, 4200, 5200);
68  m_hHit = new TH2F("hHit", "hHit", 56, 0, 56, 400, 0, 400);
69 
70  oldDir->cd();
71  }
72 
74  {
75  REG_HISTOGRAM
76  StoreArray<CDCHit> m_cdcHits;
77  StoreArray<CDCRawHit> m_cdcRawHits;
78  StoreObjPtr<TRGSummary> m_trgSummary;
79  m_cdcHits.isOptional();
80  m_cdcRawHits.isOptional();
81  m_trgSummary.isOptional();
82  }
83 
85  {
86 
87  m_hNEvents->Reset();
88  m_hADC->Reset();
89  m_hADCTOTCut->Reset();
90  m_hTDC->Reset();
91  m_hHit->Reset();
92  m_hOcc->Reset();
93  }
94 
96  {
97  StoreArray<CDCHit> m_cdcHits;
98  StoreArray<CDCRawHit> m_cdcRawHits;
99  StoreObjPtr<TRGSummary> m_trgSummary;
100  const int nWires = 14336;
101  setReturnValue(1);
102  if (!m_trgSummary.isValid() || (m_trgSummary->getTimType() == Belle2::TRGSummary::TTYP_RAND)) {
103  setReturnValue(0);
104  return;
105  }
106 
107  if (m_cdcHits.getEntries() < m_minHits) {
108 
109  setReturnValue(0); return;
110  }
111  m_nEvents += 1;
112  m_hOcc->Fill(static_cast<float>(m_cdcHits.getEntries()) / nWires);
113  for (const auto& raw : m_cdcRawHits) {
114  int bid = raw.getBoardId();
115  int adc = raw.getFADC();
116  int tdc = raw.getTDC();
117  int tot = raw.getTOT();
118  m_hADC->Fill(bid, adc);
119  if (tot > 4) {
120  m_hADCTOTCut->Fill(bid, adc);
121  }
122 
123  if (adc > 50 && tot > 1) {
124  m_hTDC->Fill(bid, tdc);
125  }
126  }
127  for (const auto& hit : m_cdcHits) {
128  int lay = hit.getICLayer();
129  int wire = hit.getIWire();
130  m_hHit->Fill(lay, wire);
131  }
132  }
133 
135  {
136  m_hNEvents->SetBinContent(1, m_nEvents);
137  }
138 
140  {
141  }
143 }
Belle2::CDCDQMModule::terminate
void terminate() override
Termination action.
Definition: CDCDQMModule.cc:139
Belle2::Module::setDescription
void setDescription(const std::string &description)
Sets the description of the module.
Definition: Module.cc:216
Belle2::CDCDQMModule::m_hADC
TH2F * m_hADC
Histogram of ADC for all boards (0-299)
Definition: CDCDQMModule.h:100
REG_MODULE
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition: Module.h:652
Belle2::Module::c_ParallelProcessingCertified
@ c_ParallelProcessingCertified
This module can be run in parallel processing mode safely (All I/O must be done through the data stor...
Definition: Module.h:82
Belle2::CDCDQMModule::m_hOcc
TH1F * m_hOcc
Histogram of occupancy.
Definition: CDCDQMModule.h:104
Belle2::TRGSummary::TTYP_RAND
@ TTYP_RAND
random trigger events
Definition: TRGSummary.h:78
Belle2::CDCDQMModule::event
void event() override
Event processor.
Definition: CDCDQMModule.cc:95
Belle2::CDCDQMModule::beginRun
void beginRun() override
Called when entering a new run.
Definition: CDCDQMModule.cc:84
Belle2::CDCDQMModule::~CDCDQMModule
virtual ~CDCDQMModule()
Destructor.
Definition: CDCDQMModule.cc:51
Belle2::Module::setPropertyFlags
void setPropertyFlags(unsigned int propertyFlags)
Sets the flags for the module properties.
Definition: Module.cc:210
Belle2::CDCDQMModule::endRun
void endRun() override
End-of-run action.
Definition: CDCDQMModule.cc:134
Belle2::CDCDQMModule::m_nEvents
Long64_t m_nEvents
Number of events processed.
Definition: CDCDQMModule.h:98
Belle2::CDCDQMModule::m_minHits
int m_minHits
Minimum hits for processing.
Definition: CDCDQMModule.h:105
Belle2::CDCDQMModule::defineHisto
void defineHisto() override
Histogram definitions.
Definition: CDCDQMModule.cc:55
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::StoreObjPtr
Type-safe access to single objects in the data store.
Definition: ParticleList.h:33
Belle2::CDCDQMModule::m_hNEvents
TH1F * m_hNEvents
Histogram of num.
Definition: CDCDQMModule.h:99
Belle2::CDCDQMModule::m_hHit
TH2F * m_hHit
Histogram of hits for all layers (0-55)
Definition: CDCDQMModule.h:103
Belle2::CDCDQMModule::m_hTDC
TH2F * m_hTDC
Histogram of TDC for all boards (0-299)
Definition: CDCDQMModule.h:102
Belle2::Module::setReturnValue
void setReturnValue(int value)
Sets the return value for this module as integer.
Definition: Module.cc:222
Belle2::Module::addParam
void addParam(const std::string &name, T &paramVariable, const std::string &description, const T &defaultValue)
Adds a new parameter to the module.
Definition: Module.h:562
Belle2::CDCDQMModule::initialize
void initialize() override
Initialize the Module.
Definition: CDCDQMModule.cc:73
Belle2::StoreArray
Accessor to arrays stored in the data store.
Definition: ECLMatchingPerformanceExpertModule.h:33
Belle2::CDCDQMModule::m_hADCTOTCut
TH2F * m_hADCTOTCut
Histogram of ADC with tot cut for all boards (0-299)
Definition: CDCDQMModule.h:101
Belle2::StoreArray::getEntries
int getEntries() const
Get the number of objects in the array.
Definition: StoreArray.h:226
Belle2::HistoModule
HistoModule.h is supposed to be used instead of Module.h for the modules with histogram definitions t...
Definition: HistoModule.h:29
Belle2::StoreObjPtr::isValid
bool isValid() const
Check whether the object was created.
Definition: StoreObjPtr.h:120