ROOT  6.07/01
Reference Guide
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Groups Pages
TEveCaloData.h
Go to the documentation of this file.
1 // @(#)root/eve:$Id$
2 // Author: Matevz Tadel 2007
3 
4 /*************************************************************************
5  * Copyright (C) 1995-2007, Rene Brun and Fons Rademakers. *
6  * All rights reserved. *
7  * *
8  * For the licensing terms see $ROOTSYS/LICENSE. *
9  * For the list of contributors see $ROOTSYS/README/CREDITS. *
10  *************************************************************************/
11 
12 #ifndef ROOT_TEveCaloData
13 #define ROOT_TEveCaloData
14 
15 #include <vector>
16 #include "TEveElement.h"
17 
18 #include "TMath.h"
19 
20 class TGLSelectRecord;
21 
22 class TH2F;
23 class TAxis;
24 class THStack;
25 
26 class TEveCaloData: public TEveElement,
27  public TNamed
28 {
29 public:
30  struct SliceInfo_t
31  {
32  TString fName; // Name of the slice, eg. ECAL, HCAL.
33  Float_t fThreshold; // Only display towers with higher energy.
34  Color_t fColor; // Color used to draw this longitudinal slice.
35  Color_t fTransparency; // Transparency used to draw this longitudinal slice.
36 
38 
39  virtual ~SliceInfo_t() {}
40 
41  void Setup(const char* name, Float_t threshold, Color_t col, Char_t transp = 101)
42  {
43  fName = name;
44  fThreshold = threshold;
45  fColor = col;
46  if (transp <= 100) fTransparency = transp;
47  };
48 
49  ClassDef(SliceInfo_t, 0); // Slice info for histogram stack.
50  };
51 
52  typedef std::vector<SliceInfo_t> vSliceInfo_t;
53  typedef std::vector<SliceInfo_t>::iterator vSliceInfo_i;
54 
55  /**************************************************************************/
56 
57  struct CellId_t
58  {
59  // Cell ID inner structure.
60 
63 
65 
67 
68  bool operator<(const CellId_t& o) const
69  { return (fTower == o.fTower) ? fSlice < o.fSlice : fTower < o.fTower; }
70  };
71 
72  struct CellGeom_t
73  {
74  // Cell geometry inner structure.
75 
80 
81  Float_t fThetaMin; // cached
82  Float_t fThetaMax; // cached
83 
85  CellGeom_t(Float_t etaMin, Float_t etaMax, Float_t phiMin, Float_t phiMax) {Configure(etaMin, etaMax, phiMin, phiMax);}
86  virtual ~CellGeom_t() {}
87 
88  void Configure(Float_t etaMin, Float_t etaMax, Float_t phiMin, Float_t phiMax);
89 
90  Float_t EtaMin() const { return fEtaMin; }
91  Float_t EtaMax() const { return fEtaMax; }
92  Float_t Eta() const { return (fEtaMin+fEtaMax)*0.5f; }
93  Float_t EtaDelta() const { return fEtaMax-fEtaMin; }
94 
95  Float_t PhiMin() const { return fPhiMin; }
96  Float_t PhiMax() const { return fPhiMax; }
97  Float_t Phi() const { return (fPhiMin+fPhiMax)*0.5f; }
98  Float_t PhiDelta() const { return fPhiMax-fPhiMin; }
99 
100  Float_t ThetaMin() const { return fThetaMin; }
101  Float_t ThetaMax() const { return fThetaMax; }
102  Float_t Theta() const { return (fThetaMax+fThetaMin)*0.5f; }
103  Float_t ThetaDelta() const { return fThetaMax-fThetaMin; }
104 
106  {
107  const Float_t phi = Phi();
108  return ((phi > 0 && phi <= TMath::Pi()) || phi < - TMath::Pi());
109  }
110 
111  virtual void Dump() const;
112  };
113 
114  struct CellData_t : public CellGeom_t
115  {
116  // Cell data inner structure.
117 
119 
121  virtual ~CellData_t() {}
122 
123  Float_t Value(Bool_t) const;
124  virtual void Dump() const;
125  };
126 
127 
128  struct RebinData_t
129  {
131 
132  std::vector<Float_t> fSliceData;
133  std::vector<Int_t> fBinData;
134 
135  Float_t* GetSliceVals(Int_t bin);
136 
137  void Clear()
138  {
139  fSliceData.clear();
140  fBinData.clear();
141  }
142  };
143 
144  /**************************************************************************/
145 
146  typedef std::vector<CellId_t> vCellId_t;
147  typedef std::vector<CellId_t>::iterator vCellId_i;
148 
149  typedef std::vector<CellGeom_t> vCellGeom_t;
150  typedef std::vector<CellGeom_t>::iterator vCellGeom_i;
151  typedef std::vector<CellGeom_t>::const_iterator vCellGeom_ci;
152 
153 private:
154  TEveCaloData(const TEveCaloData&); // Not implemented
155  TEveCaloData& operator=(const TEveCaloData&); // Not implemented
156 
157 protected:
159 
162 
164 
165  Float_t fMaxValEt; // cached
166  Float_t fMaxValE; // cached
167 
169 
172 
173 public:
174  TEveCaloData(const char* n="TEveCalData", const char* t="");
175  virtual ~TEveCaloData() {}
176 
177  virtual void UnSelected();
178  virtual void UnHighlighted();
179 
180  virtual TString GetHighlightTooltip();
181 
182  virtual void FillImpliedSelectedSet(Set_t& impSelSet);
183 
184  virtual void GetCellList(Float_t etaMin, Float_t etaMax,
185  Float_t phi, Float_t phiRng,
186  vCellId_t &out) const = 0;
187 
190  void PrintCellsSelected();
191  void ProcessSelection(vCellId_t& sel_cells, TGLSelectRecord& rec);
192 
193  virtual void Rebin(TAxis *ax, TAxis *ay, vCellId_t &in, Bool_t et, RebinData_t &out) const = 0;
194 
195 
196  virtual void GetCellData(const CellId_t &id, CellData_t& data) const = 0;
197 
198  virtual void InvalidateUsersCellIdCache();
199  virtual void DataChanged();
200  virtual void CellSelectionChanged();
201 
202  Int_t GetNSlices() const { return fSliceInfos.size(); }
204  void SetSliceThreshold(Int_t slice, Float_t threshold);
205  Float_t GetSliceThreshold(Int_t slice) const;
206  void SetSliceColor(Int_t slice, Color_t col);
207  Color_t GetSliceColor(Int_t slice) const;
208  void SetSliceTransparency(Int_t slice, Char_t t);
209  Char_t GetSliceTransparency(Int_t slice) const;
210 
211  virtual void GetEtaLimits(Double_t &min, Double_t &max) const = 0;
212 
213  virtual void GetPhiLimits(Double_t &min, Double_t &max) const = 0;
214 
215  virtual Float_t GetMaxVal(Bool_t et) const { return et ? fMaxValEt : fMaxValE; }
216  Bool_t Empty() const { return fMaxValEt < 1e-5; }
217 
218  virtual TAxis* GetEtaBins() const { return fEtaAxis; }
219  virtual void SetEtaBins(TAxis* ax) { fEtaAxis=ax; }
220 
221  virtual TAxis* GetPhiBins() const { return fPhiAxis; }
222  virtual void SetPhiBins(TAxis* ax) { fPhiAxis=ax; }
223 
224  virtual Float_t GetEps() const { return fEps; }
225  virtual void SetEps(Float_t eps) { fEps=eps; }
226 
227  Bool_t GetWrapTwoPi() const { return fWrapTwoPi; }
229 
230  static Float_t EtaToTheta(Float_t eta);
231 
232 
233  ClassDef(TEveCaloData, 0); // Manages calorimeter event data.
234 };
235 
236 /**************************************************************************/
237 /**************************************************************************/
238 
240 {
241 
242 private:
243  TEveCaloDataVec(const TEveCaloDataVec&); // Not implemented
244  TEveCaloDataVec& operator=(const TEveCaloDataVec&); // Not implemented
245 
246 protected:
247  typedef std::vector<Float_t> vFloat_t;
248  typedef std::vector<Float_t>::iterator vFloat_i;
249 
250  typedef std::vector<vFloat_t> vvFloat_t;
251  typedef std::vector<vFloat_t>::iterator vvFloat_i;
252 
255 
256  Int_t fTower; // current tower
257 
260 
263 
264 public:
265  TEveCaloDataVec(Int_t nslices);
266  virtual ~TEveCaloDataVec();
267 
268  Int_t AddSlice();
269  Int_t AddTower(Float_t etaMin, Float_t etaMax, Float_t phiMin, Float_t phiMax);
270  void FillSlice(Int_t slice, Float_t value);
271  void FillSlice(Int_t slice, Int_t tower, Float_t value);
272 
273  Int_t GetNCells() { return fGeomVec.size(); }
274  std::vector<Float_t>& GetSliceVals(Int_t slice) { return fSliceVec[slice]; }
275  std::vector<TEveCaloData::CellGeom_t>& GetCellGeom() { return fGeomVec; }
276 
277  virtual void GetCellList(Float_t etaMin, Float_t etaMax,
278  Float_t phi, Float_t phiRng,
279  vCellId_t &out) const;
280 
281  virtual void Rebin(TAxis *ax, TAxis *ay, vCellId_t &in, Bool_t et, RebinData_t &out) const;
282 
283  virtual void GetCellData(const TEveCaloData::CellId_t &id, TEveCaloData::CellData_t& data) const;
284  virtual void GetEtaLimits(Double_t &min, Double_t &max) const { min=fEtaMin, max=fEtaMax;}
285  virtual void GetPhiLimits(Double_t &min, Double_t &max) const { min=fPhiMin; max=fPhiMax;}
286 
287 
288  virtual void DataChanged();
289  void SetAxisFromBins(Double_t epsX=0.001, Double_t epsY=0.001);
290 
291  ClassDef(TEveCaloDataVec, 0); // Manages calorimeter event data.
292 };
293 
294 /**************************************************************************/
295 /**************************************************************************/
296 
298 {
299 private:
300  TEveCaloDataHist(const TEveCaloDataHist&); // Not implemented
301  TEveCaloDataHist& operator=(const TEveCaloDataHist&); // Not implemented
302 
303 protected:
305 
306 public:
308  virtual ~TEveCaloDataHist();
309 
310  virtual void GetCellList( Float_t etaMin, Float_t etaMax,
311  Float_t phi, Float_t phiRng, vCellId_t &out) const;
312 
313  virtual void Rebin(TAxis *ax, TAxis *ay, vCellId_t &in, Bool_t et, RebinData_t &out) const;
314 
315  virtual void GetCellData(const TEveCaloData::CellId_t &id, TEveCaloData::CellData_t& data) const;
316 
317  virtual void GetEtaLimits(Double_t &min, Double_t &max) const;
318  virtual void GetPhiLimits(Double_t &min, Double_t &max) const;
319 
320 
321  virtual void DataChanged();
322 
323  THStack* GetStack() { return fHStack; }
324 
325  TH2F* GetHist(Int_t slice) const;
326 
327  Int_t AddHistogram(TH2F* hist);
328 
329  ClassDef(TEveCaloDataHist, 0); // Manages calorimeter TH2F event data.
330 };
331 
332 #endif
333 
Float_t GetSliceThreshold(Int_t slice) const
Get threshold for given slice.
ClassDef(TEveCaloDataVec, 0)
Bool_t IsUpperRho() const
Definition: TEveCaloData.h:105
vCellId_t & GetCellsSelected()
Definition: TEveCaloData.h:188
TAxis * fEtaAxis
Definition: TEveCaloData.h:160
virtual TAxis * GetPhiBins() const
Definition: TEveCaloData.h:221
A central manager for calorimeter event data.
Definition: TEveCaloData.h:26
static Vc_ALWAYS_INLINE int_v min(const int_v &x, const int_v &y)
Definition: vector.h:433
std::vector< CellId_t >::iterator vCellId_i
Definition: TEveCaloData.h:147
The Histogram stack class.
Definition: THStack.h:35
SliceInfo_t & RefSliceInfo(Int_t s)
Definition: TEveCaloData.h:203
TEveCaloDataVec(const TEveCaloDataVec &)
void Setup(const char *name, Float_t threshold, Color_t col, Char_t transp=101)
Definition: TEveCaloData.h:41
vCellId_t fCellsHighlighted
Definition: TEveCaloData.h:171
virtual Float_t GetMaxVal(Bool_t et) const
Definition: TEveCaloData.h:215
float Float_t
Definition: RtypesCore.h:53
virtual Float_t GetEps() const
Definition: TEveCaloData.h:224
Color_t GetSliceColor(Int_t slice) const
Get color for given slice.
std::set< TEveElement * > Set_t
Definition: TEveElement.h:73
void Configure(Float_t etaMin, Float_t etaMax, Float_t phiMin, Float_t phiMax)
Definition: Rtypes.h:61
std::vector< Int_t > fBinData
Definition: TEveCaloData.h:133
ClassDef(TEveCaloDataHist, 0)
virtual void DataChanged()
Tell users (TEveCaloViz instances using this data) that data has changed and they should update the l...
virtual void GetEtaLimits(Double_t &min, Double_t &max) const
Get eta limits.
Char_t GetSliceTransparency(Int_t slice) const
Get transparency for given slice.
std::vector< vFloat_t > vvFloat_t
Definition: TEveCaloData.h:250
virtual ~TEveCaloData()
Definition: TEveCaloData.h:175
virtual TString GetHighlightTooltip()
Float_t PhiMax() const
Definition: TEveCaloData.h:96
std::vector< vFloat_t >::iterator vvFloat_i
Definition: TEveCaloData.h:251
Calo data for universal cell geometry.
Definition: TEveCaloData.h:239
Basic string class.
Definition: TString.h:137
virtual void GetCellData(const TEveCaloData::CellId_t &id, TEveCaloData::CellData_t &data) const
Get cell geometry and value from cell ID.
Bool_t GetWrapTwoPi() const
Definition: TEveCaloData.h:227
TEveCaloData(const TEveCaloData &)
int Int_t
Definition: RtypesCore.h:41
void FillSlice(Int_t slice, Float_t value)
Fill given slice in the current tower.
bool Bool_t
Definition: RtypesCore.h:59
T etaMax()
Function providing the maximum possible value of pseudorapidity for a non-zero rho, in the Scalar type with the largest dynamic range.
Definition: etaMax.h:50
virtual void Rebin(TAxis *ax, TAxis *ay, vCellId_t &in, Bool_t et, RebinData_t &out) const
Rebin.
void SetSliceColor(Int_t slice, Color_t col)
Set color for given slice.
Bool_t fWrapTwoPi
Definition: TEveCaloData.h:163
std::vector< Float_t > & GetSliceVals(Int_t slice)
Definition: TEveCaloData.h:274
TFile * f
virtual void UnHighlighted()
Virtual method TEveElement::UnHighlighted.
Float_t PhiMin() const
Definition: TEveCaloData.h:95
virtual void FillImpliedSelectedSet(Set_t &impSelSet)
Populate set impSelSet with derived / dependant elements.
Float_t fEps
Definition: TEveCaloData.h:168
virtual void GetCellList(Float_t etaMin, Float_t etaMax, Float_t phi, Float_t phiRng, vCellId_t &out) const =0
std::vector< SliceInfo_t >::iterator vSliceInfo_i
Definition: TEveCaloData.h:53
static Float_t EtaToTheta(Float_t eta)
Float_t PhiDelta() const
Definition: TEveCaloData.h:98
Int_t GetNSlices() const
Definition: TEveCaloData.h:202
The TNamed class is the base class for all named ROOT classes.
Definition: TNamed.h:33
TAxis * fPhiAxis
Definition: TEveCaloData.h:161
std::vector< CellGeom_t >::iterator vCellGeom_i
Definition: TEveCaloData.h:150
virtual ~TEveCaloDataHist()
Destructor.
Float_t Theta() const
Definition: TEveCaloData.h:102
virtual void InvalidateUsersCellIdCache()
Invalidate cell ids cache on back ptr references.
std::vector< SliceInfo_t > vSliceInfo_t
Definition: TEveCaloData.h:52
std::vector< CellId_t > vCellId_t
Definition: TEveCaloData.h:146
Float_t EtaMin() const
Definition: TEveCaloData.h:90
virtual void GetPhiLimits(Double_t &min, Double_t &max) const =0
std::vector< Float_t > vFloat_t
Definition: TEveCaloData.h:247
char * out
Definition: TBase64.cxx:29
short Color_t
Definition: RtypesCore.h:79
virtual void UnSelected()
Virtual method TEveElement::UnSelect.
Float_t ThetaMin() const
Definition: TEveCaloData.h:100
Float_t fMaxValE
Definition: TEveCaloData.h:166
Int_t AddHistogram(TH2F *hist)
Add new slice to calo tower.
virtual void SetEps(Float_t eps)
Definition: TEveCaloData.h:225
Float_t ThetaMax() const
Definition: TEveCaloData.h:101
CellId_t(Int_t t, Int_t s, Float_t f=1.0f)
Definition: TEveCaloData.h:66
TThread * t[5]
Definition: threadsh1.C:13
std::vector< CellGeom_t > vCellGeom_t
Definition: TEveCaloData.h:149
vCellId_t fCellsSelected
Definition: TEveCaloData.h:170
Class to manage histogram axis.
Definition: TAxis.h:36
virtual void DataChanged()
Update limits and notify data users.
void ProcessSelection(vCellId_t &sel_cells, TGLSelectRecord &rec)
Process newly selected cells with given select-record.
TEveCaloDataHist()
Constructor.
virtual void GetCellData(const TEveCaloData::CellId_t &id, TEveCaloData::CellData_t &data) const
Get cell geometry and value from cell ID.
std::vector< Float_t >::iterator vFloat_i
Definition: TEveCaloData.h:248
2-D histogram with a float per channel (see TH1 documentation)}
Definition: TH2.h:256
Standard selection record including information about containing scene and details ob out selected ob...
bool operator<(const CellId_t &o) const
Definition: TEveCaloData.h:68
Bool_t Empty() const
Definition: TEveCaloData.h:216
tuple w
Definition: qtexample.py:51
Cell geometry inner structure.
Definition: TEveCaloData.h:72
virtual void GetCellData(const CellId_t &id, CellData_t &data) const =0
virtual void CellSelectionChanged()
Tell users (TEveCaloViz instances using this data) that cell selection has changed and they should up...
virtual void Dump() const
Print member data.
CellGeom_t(Float_t etaMin, Float_t etaMax, Float_t phiMin, Float_t phiMax)
Definition: TEveCaloData.h:85
TEveCaloDataHist & operator=(const TEveCaloDataHist &)
ClassDef(SliceInfo_t, 0)
void SetAxisFromBins(Double_t epsX=0.001, Double_t epsY=0.001)
Set XY axis from cells geometry.
vCellId_t & GetCellsHighlighted()
Definition: TEveCaloData.h:189
std::vector< TEveCaloData::CellGeom_t > & GetCellGeom()
Definition: TEveCaloData.h:275
virtual void GetEtaLimits(Double_t &min, Double_t &max) const
Definition: TEveCaloData.h:284
Float_t Value(Bool_t) const
Return energy value associated with the cell, usually Et.
virtual void GetCellList(Float_t etaMin, Float_t etaMax, Float_t phi, Float_t phiRng, vCellId_t &out) const
Get list of cell IDs in given eta and phi range.
Float_t fMaxValEt
Definition: TEveCaloData.h:165
TEveCaloData & operator=(const TEveCaloData &)
Double_t Pi()
Definition: TMath.h:44
Float_t phi
Definition: shapesAnim.C:6
Int_t AddTower(Float_t etaMin, Float_t etaMax, Float_t phiMin, Float_t phiMax)
Add tower within eta/phi range.
Float_t EtaMax() const
Definition: TEveCaloData.h:91
vSliceInfo_t fSliceInfos
Definition: TEveCaloData.h:158
ClassDef(TEveCaloData, 0)
virtual void GetCellList(Float_t etaMin, Float_t etaMax, Float_t phi, Float_t phiRng, vCellId_t &out) const
Get list of cell-ids for given eta/phi range.
double Double_t
Definition: RtypesCore.h:55
void SetWrapTwoPi(Bool_t w)
Definition: TEveCaloData.h:228
virtual void DataChanged()
Update limits and notify data users.
vvFloat_t fSliceVec
Definition: TEveCaloData.h:253
void SetSliceTransparency(Int_t slice, Char_t t)
Set transparency for given slice.
TH2F * GetHist(Int_t slice) const
Get histogram in given slice.
static Vc_ALWAYS_INLINE int_v max(const int_v &x, const int_v &y)
Definition: vector.h:440
#define name(a, b)
Definition: linkTestLib0.cpp:5
virtual void Rebin(TAxis *ax, TAxis *ay, vCellId_t &in, Bool_t et, RebinData_t &out) const
Rebin cells.
std::vector< Float_t > fSliceData
Definition: TEveCaloData.h:132
std::vector< CellGeom_t >::const_iterator vCellGeom_ci
Definition: TEveCaloData.h:151
char Char_t
Definition: RtypesCore.h:29
virtual void GetEtaLimits(Double_t &min, Double_t &max) const =0
Float_t ThetaDelta() const
Definition: TEveCaloData.h:103
Float_t EtaDelta() const
Definition: TEveCaloData.h:93
virtual TAxis * GetEtaBins() const
Definition: TEveCaloData.h:218
virtual void GetPhiLimits(Double_t &min, Double_t &max) const
Get phi limits.
virtual ~TEveCaloDataVec()
Destructor.
TEveCaloDataVec & operator=(const TEveCaloDataVec &)
virtual void GetPhiLimits(Double_t &min, Double_t &max) const
Definition: TEveCaloData.h:285
Float_t Phi() const
Definition: TEveCaloData.h:97
Cell data inner structure.
Definition: TEveCaloData.h:114
vCellGeom_t fGeomVec
Definition: TEveCaloData.h:254
void PrintCellsSelected()
Print selected cells info.
virtual void Dump() const
Print member data.
THStack * fHStack
Definition: TEveCaloData.h:304
THStack * GetStack()
Definition: TEveCaloData.h:323
void SetSliceThreshold(Int_t slice, Float_t threshold)
Set threshold for given slice.
float value
Definition: math.cpp:443
Base class for TEveUtil visualization elements, providing hierarchy management, rendering control and...
Definition: TEveElement.h:33
virtual void SetEtaBins(TAxis *ax)
Definition: TEveCaloData.h:219
virtual void Rebin(TAxis *ax, TAxis *ay, vCellId_t &in, Bool_t et, RebinData_t &out) const =0
const Int_t n
Definition: legend1.C:16
Float_t Eta() const
Definition: TEveCaloData.h:92
virtual void SetPhiBins(TAxis *ax)
Definition: TEveCaloData.h:222
Float_t * GetSliceVals(Int_t bin)
void transp()
Definition: transp.C:17
A central manager for calorimeter data of an event written in TH2F.
Definition: TEveCaloData.h:297
Int_t AddSlice()
Add new slice.