Belle II Software  release-05-01-25
TRGGRL.cc
1 //-----------------------------------------------------------------------------
2 // $Id$
3 //-----------------------------------------------------------------------------
4 // Filename : TRGGRL.cc
5 // Section : TRG GRL
6 // Owner : Yun-Tsung Lai
7 // Email : ytlai@hep1.phys.ntu.edu.tw
8 //-----------------------------------------------------------------------------
9 // Description : A class to represent GRL.
10 //-----------------------------------------------------------------------------
11 // 0.00 : 2013/12/13 : First version
12 //-----------------------------------------------------------------------------
13 
14 #define TRG_SHORT_NAMES
15 #define TRGGRL_SHORT_NAMES
16 
17 #include <fstream>
18 #include "framework/datastore/StoreArray.h"
19 #include "trg/trg/Utilities.h"
20 #include "trg/grl/TRGGRL.h"
21 #include "trg/cdc/TRGCDC.h"
22 #include "trg/cdc/Track.h"
23 #include "trg/ecl/dataobjects/TRGECLCluster.h"
24 #include "trg/ecl/TrgEclCluster.h"
25 #include <math.h>
26 #include <TFile.h>
27 #include <TTree.h>
28 #include <framework/logging/Logger.h>
29 
30 # define M_PI 3.14159265358979323846
31 
32 using namespace std;
33 
34 namespace Belle2 {
40  TRGGRL*
41  TRGGRL::_grl = 0;
42 
43  string
44  TRGGRL::name(void) const
45  {
46  return "TRGGRL";
47  }
48 
49  string
50  TRGGRL::version(void) const
51  {
52  return string("TRGGRL 0.01");
53  }
54 
55  TRGGRL*
56  TRGGRL::getTRGGRL(const string& configFile,
57  unsigned simulationMode,
58  unsigned fastSimulationMode,
59  unsigned firmwareSimulationMode)
60  {
61  if (_grl) {
62  //delete _grl;
63  _grl = 0;
64  }
65 
66  if (configFile != "good-bye") {
67  _grl = new TRGGRL(configFile,
68  simulationMode,
69  fastSimulationMode,
70  firmwareSimulationMode);
71  } else {
72  B2DEBUG(100, "TRGGRL::getTRGGRL ... good-bye");
73  _grl = 0;
74  }
75 
76  return _grl;
77  }
78 
79  TRGGRL*
80  TRGGRL::getTRGGRL(void)
81  {
82  if (! _grl)
83  B2WARNING("TRGGRL::getTRGGRL !!! TRGGRL is not created yet");
84  return _grl;
85  }
86 
87  TRGGRL::TRGGRL(const string& configFile,
88  unsigned simulationMode,
89  unsigned fastSimulationMode,
90  unsigned firmwareSimulationMode)
91  : _debugLevel(0),
92  _configFilename(configFile),
93  _simulationMode(simulationMode),
94  _fastSimulationMode(fastSimulationMode),
95  _firmwareSimulationMode(firmwareSimulationMode),
96  _clock(Belle2_GDL::GDLSystemClock)
97  {
98 
99  B2DEBUG(100, "TRGGRL ... TRGGRL initializing with " << _configFilename
100  << " mode=0x" << hex << _simulationMode << dec);
101 
102  initialize();
103 
104  B2DEBUG(100, "TRGGRL ... TRGGRL created with " << _configFilename);
105  }
106 
107  void
109  {
110 
111  m_file = new TFile("trggrl.root", "RECREATE");
112  h1 = new TTree("h1", "h1");
113 
114  h1->Branch("3d", &x0, "3d");
115  h1->Branch("dr", &x1, "dr");
116  h1->Branch("dz", &x2, "dz");
117  h1->Branch("poe", &x3, "poe");
118  h1->Branch("z0", &x4, "z0");
119  h1->Branch("pt", &x5, "pt");
120  h1->Branch("pz", &x6, "pz");
121  h1->Branch("e", &x7, "e");
122 
123  configure();
124  }
125 
126  void
128  {
129  }
130 
131  void
132  TRGGRL::dump(const string& msg) const
133  {
134 
135  if (msg != "") B2DEBUG(100, "dump nothing...");
136 
137  }
138 
139  void
141  {
142  }
143 
144  void
146  {
147  }
148 
149  void
151  {
152 
153  matchList.clear();
154  matchList3D.clear();
155 
156  B2DEBUG(100, "do nothing...");
157 
158  }
159 
161  {
162  h1->Write();
163  m_file->Write();
164  m_file->Close();
165  clear();
166  }
167 
168  void
170  {
171 
173  vector<TRGCDCTrack*> trackList = _cdc->getTrackList2D();
174  vector<TRGCDCTrack*> trackList3D = _cdc->getTrackList3D();
175  StoreArray<TRGECLCluster> ClusterArray;
176 
177  unsigned n_track = trackList.size();
178  unsigned n_track3D = trackList3D.size();
179  unsigned n_cluster = ClusterArray.getEntries();
180 
181 // if (TRGDebug::level() > 2) cout <<"yt_grl "<< n_cluster << " " << n_track << endl;
182 
183  for (unsigned i = 0; i < n_track; i++) {
184  // vector<TRGGRLMatch *> match_i;
185  if (n_cluster == 0) break;
186  else if (n_cluster == 1) {
187  TRGGRLMatch* match = new TRGGRLMatch(trackList[i], ClusterArray[0], 0);
188  matchList.push_back(match);
189  } else if (n_cluster > 1) {
190  int best_j = 0; double old_dr = 99999;
191  for (unsigned j = 0; j < n_cluster; j++) {
192  TRGGRLMatch* match = new TRGGRLMatch(trackList[i], ClusterArray[j], 0);
193  if (match->getDr() < old_dr) {best_j = j; old_dr = match->getDr();}
194  }
195  TRGGRLMatch* match = new TRGGRLMatch(trackList[i], ClusterArray[best_j], 0);
196  matchList.push_back(match);
197  }
198  }
199 
200  for (unsigned i = 0; i < n_track3D; i++) {
201  // vector<TRGGRLMatch *> match_i;
202  if (n_cluster == 0) break;
203  else if (n_cluster == 1) {
204  TRGGRLMatch* match = new TRGGRLMatch(trackList3D[i], ClusterArray[0], 1);
205  matchList.push_back(match);
206  } else if (n_cluster > 1) {
207  int best_j = 0; double old_dr = 99999;
208  for (unsigned j = 0; j < n_cluster; j++) {
209  TRGGRLMatch* match = new TRGGRLMatch(trackList3D[i], ClusterArray[j], 0);
210  if (match->getDr() < old_dr) {best_j = j; old_dr = match->getDr();}
211  }
212  TRGGRLMatch* match = new TRGGRLMatch(trackList3D[i], ClusterArray[best_j], 0);
213  matchList3D.push_back(match);
214  }
215  }
216 
217  //-----Fill tree
218  for (unsigned i = 0; i < matchList.size(); i++) {
219  TRGGRLMatch& match = * matchList[i];
220 
221  x0 = match.getMatch3D();
222  x1 = match.getDr();
223  x2 = match.getDz();
224  x3 = match.getPoe();
225  x4 = match.getCenter_z0();
226  x5 = match.getCenter_pt();
227  x6 = match.getCenter_pz();
228  x7 = match.getCluster_e();
229  h1->Fill();
230  }
231 
232  for (unsigned i = 0; i < matchList3D.size(); i++) {
233  TRGGRLMatch& match = * matchList3D[i];
234 
235  x0 = match.getMatch3D();
236  x1 = match.getDr();
237  x2 = match.getDz();
238  x3 = match.getPoe();
239  x4 = match.getCenter_z0();
240  x5 = match.getCenter_pt();
241  x6 = match.getCenter_pz();
242  x7 = match.getCluster_e();
243  h1->Fill();
244  }
245 
246  //--------------------------------------
247 
248 
249  const bool fast = (_simulationMode & 1);
250  const bool firm = (_simulationMode & 2);
251  if (fast)
252  fastSimulation();
253  if (firm)
255  }
256 
257  void
259  {
260  }
261 
262  void
264  {
265  }
266 
267  void
269  {
270  }
271  /*
272  bool
273  TRGGRL::barrel_matching_2D(TRGCDCTrack * track, TRGECLCluster * cluster){
274 
275  //-- track/TRGCDC information
276  const TRGCDCHelix & helix = track->helix();
277  double pt = track->pt();
278  double center_x = helix.center().x(), center_y = helix.center().y(), center_z = helix.center().z();
279  double r = sqrt(center_x*center_x + center_y*center_y);//helix.radius();
280  double phi = atan2(center_y, center_x) ;
281 
282  //-- cluster/TRGECL information
283  double cluster_x = cluster->getPositionX(), cluster_y = cluster->getPositionY(), cluster_z = cluster->getPositionZ();
284  double cluster_e = cluster->getEnergyDep();
285  double R = sqrt(cluster_x*cluster_x + cluster_y*cluster_y);
286  double D = sqrt(cluster_x*cluster_x + cluster_y*cluster_y + cluster_z*cluster_z);
287  double re_scaled_p = pt*D/R;
288 
289  //-- calculation
290 
291  double theta0 = acos(R/(2*r)) + phi;
292  double theta1 = 2*phi - theta0;
293 
294  double ex_x0 = R*cos(theta0), ex_y0 = R*sin(theta0), ex_x1 = R*cos(theta1), ex_y1 = R*sin(theta1);
295 
296  double dr0 = sqrt( (ex_x0-cluster_x)*(ex_x0-cluster_x) + (ex_y0-cluster_y)*(ex_y0-cluster_y) );
297  double dr1 = sqrt( (ex_x1-cluster_x)*(ex_x1-cluster_x) + (ex_y1-cluster_y)*(ex_y1-cluster_y) );
298  double dr = (dr0 < dr1) ? dr0 : dr1;
299 
300  if (TRGDebug::level() > 1) printf("%s %f %f %f %f %f %f %f \n","dump! ",dr,dr0,dr1,pt,re_scaled_p,cluster_e,re_scaled_p/cluster_e);
301 
302  if (TRGDebug::level() > 1 && dr > 30) {
303  cout << " " << endl;
304  cout << "double center_x = " << center_x << ";" <<endl;
305  cout << "double center_y = " << center_y << ";" <<endl;
306  cout << "double center_z = " << center_z << ";" <<endl;
307  cout << "double radius = " << r << ";" <<endl;
308  cout << "double pt = " << pt << ";" <<endl;
309  cout << "double cluster_x = " << cluster_x << ";" <<endl;
310  cout << "double cluster_y = " << cluster_y << ";" <<endl;
311  cout << "double ex_x0 = " << ex_x0 << ";" <<endl;
312  cout << "double ex_y0 = " << ex_y0 << ";" <<endl;
313  cout << "double ex_x1 = " << ex_x1 << ";" <<endl;
314  cout << "double ex_y1 = " << ex_y1 << ";" <<endl;
315  }
316 
317  return true;
318  }
319  */
320 
321 
322 
324 } // namespace Belle2
Belle2::TRGGRL
a class for TRGGRL
Definition: TRGGRL.h:48
Belle2::TRGGRL::h1
TTree * h1
root tree
Definition: TRGGRL.h:170
Belle2::TRGCDC::getTrackList3D
std::vector< TRGCDCTrack * > getTrackList3D(void)
returns 3D track list (fitted).
Definition: TRGCDC.cc:210
Belle2::TRGGRL::matchList
std::vector< TRGGRLMatch * > matchList
Vector which stores list of TRGGRLMatch without 3D information.
Definition: TRGGRL.h:193
Belle2::TRGGRL::x6
double x6
Temporary variables to make tree in root files.
Definition: TRGGRL.h:185
Belle2::TRGGRL::firmwareSimulation
void firmwareSimulation(void)
Firmware simulation.
Definition: TRGGRL.cc:263
Belle2::TRGGRL::x5
double x5
Temporary variables to make tree in root files.
Definition: TRGGRL.h:183
Belle2::TRGGRL::dump
void dump(const std::string &message) const
dumps debug information.
Definition: TRGGRL.cc:132
Belle2::TRGGRL::simulate
void simulate(void)
fast trigger simulation.
Definition: TRGGRL.cc:169
Belle2::TRGGRL::fastClear
void fastClear(void)
clears TRGGRL information.
Definition: TRGGRL.cc:145
Belle2::TRGGRL::m_file
TFile * m_file
root file
Definition: TRGGRL.h:167
Belle2::TRGGRL::x3
double x3
Temporary variables to make tree in root files.
Definition: TRGGRL.h:179
Belle2::TRGGRL::x2
double x2
Temporary variables to make tree in root files.
Definition: TRGGRL.h:177
Belle2::TRGGRL::x4
double x4
Temporary variables to make tree in root files.
Definition: TRGGRL.h:181
Belle2::TRGCDC::getTrackList2D
std::vector< TRGCDCTrack * > getTrackList2D(void)
returns 2D track list (no fit).
Definition: TRGCDC.cc:198
Belle2::TRGGRL::update
void update(bool mcAnalysis=true)
updates TRGGRL information.
Definition: TRGGRL.cc:150
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::TRGGRL::_configFilename
std::string _configFilename
root file name.
Definition: TRGGRL.h:152
Belle2::TRGCDC::getTRGCDC
static TRGCDC * getTRGCDC(void)
returns TRGCDC object.
Definition: TRGCDC.cc:190
Belle2::TRGGRL::x1
double x1
Temporary variables to make tree in root files.
Definition: TRGGRL.h:175
Belle2::TRGGRL::~TRGGRL
virtual ~TRGGRL()
Destructor.
Definition: TRGGRL.cc:160
Belle2::TRGGRL::matchList3D
std::vector< TRGGRLMatch * > matchList3D
Vector which stores list of TRGGRLMatch with 3D information.
Definition: TRGGRL.h:195
Belle2::TRGCDC
The instance of TRGCDC is a singleton.
Definition: TRGCDC.h:70
Belle2::TRGGRL::x0
double x0
Temporary variables to make tree in root files.
Definition: TRGGRL.h:173
Belle2::TRGGRLMatch
A class to represent a matching candidate in TRGGRL A matching candidate consists of a TRGCDCTrack an...
Definition: TRGGRLMatch.h:27
Belle2::TRGGRL::clear
void clear(void)
clears all TRGGRL information.
Definition: TRGGRL.cc:140
Belle2::TRGGRL::_simulationMode
unsigned _simulationMode
Simulation mode.
Definition: TRGGRL.h:155
Belle2::TRGGRL::fastSimulation
void fastSimulation(void)
Fast simulation.
Definition: TRGGRL.cc:258
Belle2::TRGGRL::x7
double x7
Temporary variables to make tree in root files.
Definition: TRGGRL.h:187
Belle2::StoreArray
Accessor to arrays stored in the data store.
Definition: ECLMatchingPerformanceExpertModule.h:33
Belle2::TRGGRL::terminate
void terminate(void)
terminates when run is finished
Definition: TRGGRL.cc:127
Belle2::TRGGRL::configure
void configure(void)
configures trigger modules for firmware simulation.
Definition: TRGGRL.cc:268
Belle2::StoreArray::getEntries
int getEntries() const
Get the number of objects in the array.
Definition: StoreArray.h:226
Belle2::TRGGRL::initialize
void initialize(void)
initializes GRL.
Definition: TRGGRL.cc:108