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