Belle II Software  release-05-02-19
ScanCDCGeoModule.cc
1 /******************************************************************
2 * BASF2 (Belle II Analysis Software Framework) *
3 * Copyright(C) YEAR Belle II Collaboration *
4 * *
5 * Test module to get CDC geometry parameters into a .root files *
6 * To run this module use script: cdc/examples/runScanGeo.py *
7 *******************************************************************/
8 #include <cdc/modules/cdcGeoScan/ScanCDCGeoModule.h>
9 #include <geometry/GeometryManager.h>
10 #include <framework/gearbox/GearDir.h>
11 #include <cdc/geometry/CDCGeometryPar.h>
12 #include "TH1F.h"
13 #include "TCanvas.h"
14 #include "TVectorF.h"
15 
16 using namespace std;
17 using namespace Belle2;
18 using namespace CDC;
19 
20 REG_MODULE(ScanCDCGeo);
21 
22 ScanCDCGeoModule::ScanCDCGeoModule()
23 {
24  setDescription("This module fills CDC geometry information in histo/tree format to a root file");
25 }
26 
27 void ScanCDCGeoModule::initialize()
28 {
29 
30  bookOutput();
31 
32  TH1F* h_nwires = new TH1F("h_nwires", "Number of wires in layer;CDC layer#", 56, 0.5, 56.5);
33  TH1F* h_stereo = new TH1F("h_stereo", "Stereo angle in this layer;CDC layer#", 56, 0.5, 56.5);
34  TH1F* h_swire_posF_phi = new TH1F("h_swire_posF_phi", "#phi of sense wires in forward pos;CDC layer#", 56, 0.5, 56.5);
35  TH1F* h_swire_posF_theta = new TH1F("h_swire_posF_theta", "#theta of sense wires in forward pos;CDC layer#", 56, 0.5, 56.5);
36  TH1F* h_swire_posB_phi = new TH1F("h_swire_posB_phi", "#phi of sense wires in backward pos;CDC layer#", 56, 0.5, 56.5);
37  TH1F* h_swire_posB_theta = new TH1F("h_swire_posB_theta", "#theta of sense wires in backward pos;CDC layer#", 56, 0.5, 56.5);
38  TH1F* h_fwire_iradius = new TH1F("h_fwire_iradius", "inner radius of field wires;CDC layer#", 56, 0.5, 56.5);
39  TH1F* h_fwire_oradius = new TH1F("h_fwire_oradius", "outer radius of field wires;CDC layer#", 56, 0.5, 56.5);
40  TH1F* h_width = new TH1F("h_width", "Cell Width;CDC layer#", 56, 0.5, 56.5);
41  TH1F* h_height = new TH1F("h_height", "Cell Height;CDC layer#", 56, 0.5, 56.5);
42  TH1F* h_length = new TH1F("h_length", "length of the wires;CDC layer#", 56, 0.5, 56.5);
43 
44  B2INFO("Creating CDCGeometryPar object");
45  CDCGeometryPar::Instance();
46  CDCGeometryPar& cdcgeo = CDCGeometryPar::Instance();
47 
48  // Print some info
49  cout << "------| Summary-1 " << endl;
50  cout << left << "Number of CDC wire layers" << setw(15) << " " << cdcgeo.nWireLayers() << endl;
51  cout << left << "Sense Wire Diameter" << setw(15) << " " << cdcgeo.senseWireDiameter() << endl;
52  cout << left << "Field Wire Diameter" << setw(15) << " " << cdcgeo.fieldWireDiameter() << endl;
53 
54  cout << "------| Summary-2 " << endl;
55  cout << left
56  << setw(4) << ""
57  << setw(10) << "Layer"
58  << setw(15) << "nwire"
59  << setw(15) << "radii sense"
60  << setw(15) << "dRdown"
61  << setw(15) << "dRup"
62  << setw(15) << "Delta"
63  << endl;
64 
65  for (unsigned int i = 0; i < cdcgeo.nWireLayers(); ++i) {
66 
67  fnWires = cdcgeo.nWiresInLayer(i);
68  h_nwires->SetBinContent(i + 1, fnWires);
69 
70  // forward position of the input sense wire
71  const TVector3& wirePosF = cdcgeo.wireForwardPosition(i, 0);
72  fswire_posF_phi = wirePosF.Phi();
73  h_swire_posF_phi->SetBinContent(i + 1, fswire_posF_phi);
74  fswire_posF_theta = wirePosF.Theta();
75  h_swire_posF_theta->SetBinContent(i + 1, fswire_posF_theta);
76 
77  const TVector3& wirePosB = cdcgeo.wireBackwardPosition(i, 0);
78  fswire_posB_phi = wirePosB.Phi();
79  h_swire_posB_phi->SetBinContent(i + 1, fswire_posB_phi);
80  fswire_posB_theta = wirePosB.Theta();
81  h_swire_posB_theta->SetBinContent(i + 1, fswire_posB_theta);
82 
83  const TVector3 wireDir = (wirePosF - wirePosB);
84  fstereoAng = wireDir.Theta();
85  if (wirePosF.Phi() < wirePosB.Phi())fstereoAng *= -1;
86  h_stereo->SetBinContent(i + 1, fstereoAng);
87 
88  fclength = wirePosF.Perp();
89  h_length->SetBinContent(i + 1, fclength);
90 
91  // the width of the cell (trapezoidal)
92  fcwidth = 2 * PI * wirePosF.Perp() / fnWires;
93  h_width->SetBinContent(i + 1, fcwidth);
94 
95  fwire_iradius = cdcgeo.innerRadiusWireLayer()[i];
96  h_fwire_iradius->SetBinContent(i + 1, fwire_iradius);
97 
98  fwire_oradius = cdcgeo.outerRadiusWireLayer()[i];
99  h_fwire_oradius->SetBinContent(i + 1, fwire_oradius);
100 
101  fcheight = fwire_oradius - fwire_iradius;
102  h_height->SetBinContent(i + 1, fcheight);
103 
104  //Printing some useful information
105  double inradiusnext = cdcgeo.innerRadiusWireLayer()[i + 1];
106  double delta = -99.0;
107  if (i < cdcgeo.nWireLayers() - 1) {
108  delta = fwire_oradius - inradiusnext;
109  }
110 
111  cout << left
112  << setw(4) << ""
113  << setw(10) << i
114  << setw(15) << fnWires
115  << setw(15) << cdcgeo.senseWireR(i)
116  << setw(15) << cdcgeo.senseWireR(i) - cdcgeo.innerRadiusWireLayer()[i]
117  << setw(15) << cdcgeo.outerRadiusWireLayer()[i] - cdcgeo.senseWireR(i)
118  << setw(15) << delta
119  << endl;
120  m_tree->Fill();
121 
122  }
123 
124  TDirectory* dhistos = m_file->mkdir("histos");
125  dhistos->cd();
126  h_nwires->Write();
127  h_swire_posF_phi->Write();
128  h_swire_posB_phi->Write();
129  h_swire_posF_theta->Write();
130  h_swire_posB_theta->Write();
131  h_stereo->Write();
132  h_fwire_iradius->Write();
133  h_fwire_oradius->Write();
134  h_width->Write();
135  h_height->Write();
136  h_length->Write();
137 
138  TDirectory* dvar = m_file->mkdir("vars");
139  dvar->cd();
140 
141  TVectorF sWireDia(1);
142  sWireDia[0] = cdcgeo.senseWireDiameter();
143  sWireDia.Write("sWireDia");
144 
145  TVectorF fWireDia(1);
146  fWireDia[0] = cdcgeo.fieldWireDiameter();
147  fWireDia.Write("fWireDia");
148 
149 }
150 
151 void ScanCDCGeoModule::bookOutput()
152 {
153  // register output root file
154  m_file = new TFile("CDCGeometryScan.root", "RECREATE");
155  m_tree = new TTree("tree", "CDC Geometry details");
156  m_tree->SetDirectory(0);
157  m_tree->Branch("lnwires", &fnWires, "lnwires/I");
158  m_tree->Branch("lsteang", &fstereoAng, "lsteang/D");
159  m_tree->Branch("lswire_fpos_phi", &fswire_posF_phi, "lswire_fpos_phi/D");
160  m_tree->Branch("lswire_fpos_theta", &fswire_posF_theta, "lswire_fpos_theta/D");
161  m_tree->Branch("lswire_bpos_phi", &fswire_posB_phi, "lswire_bpos_phi/D");
162  m_tree->Branch("lswire_bpos_theta", &fswire_posB_theta, "lswire_bpos_theta/D");
163  m_tree->Branch("lfwire_inr", &fwire_iradius, "lfwire_inr/D");
164  m_tree->Branch("lfwire_or", &fwire_oradius, "lfwire_or/D");
165  m_tree->Branch("lcwidth", &fcwidth, "lcwidth/D");
166  m_tree->Branch("lclength", &fclength, "lclength/D");
167  m_tree->Branch("lcheight", &fcheight, "lcheight/D");
168 }
169 
170 void ScanCDCGeoModule::terminate()
171 {
172  m_file->cd();
173  m_tree->Write();
174  m_file->Close();
175 }
Belle2::CDC::CDCGeometryPar::nWireLayers
unsigned nWireLayers() const
Returns a number of wire layers.
Definition: CDCGeometryPar.h:1311
Belle2::CDC::CDCGeometryPar::wireForwardPosition
const TVector3 wireForwardPosition(int layerId, int cellId, EWirePosition set=c_Base) const
Returns the forward position of the input sense wire.
Definition: CDCGeometryPar.cc:1625
Belle2::CDC::CDCGeometryPar::innerRadiusWireLayer
const double * innerRadiusWireLayer() const
Returns an array of inner radius of wire layers.
Definition: CDCGeometryPar.cc:1709
REG_MODULE
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition: Module.h:652
Belle2::CDC::CDCGeometryPar::fieldWireDiameter
double fieldWireDiameter() const
Returns diameter of the field wire.
Definition: CDCGeometryPar.h:1306
Belle2::CDC::CDCGeometryPar
The Class for CDC Geometry Parameters.
Definition: CDCGeometryPar.h:75
Belle2::CDC::CDCGeometryPar::senseWireR
double senseWireR(int layerId) const
Returns radius of sense wire in each layer.
Definition: CDCGeometryPar.h:1231
Belle2::CDC::CDCGeometryPar::wireBackwardPosition
const TVector3 wireBackwardPosition(int layerId, int cellId, EWirePosition set=c_Base) const
Returns the backward position of the input sense wire.
Definition: CDCGeometryPar.cc:1662
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::CDC::CDCGeometryPar::outerRadiusWireLayer
const double * outerRadiusWireLayer() const
Returns an array of outer radius of wire layers.
Definition: CDCGeometryPar.cc:1721
Belle2::CDC::CDCGeometryPar::nWiresInLayer
unsigned nWiresInLayer(int layerId) const
Returns wire numbers in a layer.
Definition: CDCGeometryPar.h:1211
Belle2::CDC::CDCGeometryPar::senseWireDiameter
double senseWireDiameter() const
Returns diameter of the sense wire.
Definition: CDCGeometryPar.h:1301