KaliVeda
Toolkit for HIC analysis
KVDetector.h
1 /***************************************************************************
2  kvdetector.h - description
3  -------------------
4  begin : Thu May 16 2002
5  copyright : (C) 2002 by J.D. Frankland
6  email : frankland@ganil.fr
7 
8 $Id: KVDetector.h,v 1.71 2009/05/22 14:45:40 ebonnet Exp $
9  ***************************************************************************/
10 
11 /***************************************************************************
12  * *
13  * This program is free software; you can redistribute it and/or modify *
14  * it under the terms of the GNU General Public License as published by *
15  * the Free Software Foundation; either version 2 of the License, or *
16  * (at your option) any later version. *
17  * *
18  ***************************************************************************/
19 
20 #ifndef KVDETECTOR_H
21 #define KVDETECTOR_H
22 
23 #ifndef KVD_RECPRC_CNXN
24 #define KVD_RECPRC_CNXN 1
25 #endif
26 #ifndef KVD_NORECPRC_CNXN
27 #define KVD_NORECPRC_CNXN 0
28 #endif
29 
30 #include "KVMaterial.h"
31 #include "KVPosition.h"
32 #include "KVList.h"
33 #include "KVNucleus.h"
34 #include "KVGeoDetectorNode.h"
35 #include "KVUniqueNameList.h"
36 #include "KVDetectorSignal.h"
37 #include "KVEvent.h"
38 
39 class KVGeoStrucElement;
40 class KVGroup;
41 class KVCalibrator;
42 class TGeoVolume;
43 class TTree;
44 class TGraph;
45 
157 class KVDetector: public KVMaterial, public KVPosition {
158 
159 private:
160  KVPosition fEWPosition;
161  KVUniqueNameList fParentStrucList;
162  KVGeoDetectorNode fNode;
163  static Int_t fDetCounter;
164  KVMaterial* fActiveLayer;
165  TString fNameOfArray;
166 
167  enum {
168  kIsAnalysed = BIT(14), //for reconstruction of particles
169  kActiveSet = BIT(15), //internal - flag set true when SetActiveLayer called
170  kUnidentifiedParticle = BIT(16), //set if detector is in an unidentified particle's list
171  kIdentifiedParticle = BIT(17), //set if detector is in an identified particle's list
172  };
173 
174  Int_t fIdentP;
175  Int_t fUnidentP;
176 
180  void SetMatrix(const TGeoHMatrix* m) override
181  {
183  }
184  void SetShape(TGeoBBox* s) override
185  {
187  }
188  TGeoHMatrix* GetMatrix() const override
189  {
190  return KVPosition::GetMatrix();
191  }
192  TGeoBBox* GetShape() const override
193  {
194  return KVPosition::GetShape();
195  }
196  TVector3 GetRandomPointOnSurface() const override
197  {
199  }
200  TVector3 GetSurfaceCentre() const override
201  {
203  }
204  TVector3 GetVolumeCentre() const override
205  {
207  }
208  TVector3 GetSurfaceNormal() const override
209  {
211  }
212  Double_t GetSurfaceArea(int npoints = 100000) const override
213  {
214  return KVPosition::GetSurfaceArea(npoints);
215  }
216  Double_t GetMisalignmentAngle() const override
217  {
219  }
220 
221  KVUniqueNameList fDetSignals;
222 
223  void remove_signal_for_calibrator(KVCalibrator* K);
224 
225  template<typename AbsorberStack>
226  Double_t get_corrected_energy(AbsorberStack* stack, KVNucleus* nuc, Double_t e, Bool_t transmission)
227  {
228  auto z = nuc->GetZ();
229  auto a = nuc->GetA();
230  enum SolType solution = kEmax;
231  if (!transmission) solution = kEmin;
232  Double_t EINC, ERES = GetEResAfterDetector();
233  if (transmission && ERES > 0.) {
236  EINC = stack->GetIncidentEnergyFromERes(z, a, ERES);
238  if (EINC < stack->GetEIncOfMaxDeltaE(z, a)) {
239  solution = kEmin;
241  }
245  }
246  EINC = stack->GetIncidentEnergy(z, a, e, solution);
248  if (EINC < 0) {
249  SetEResAfterDetector(-1.);
251  return EINC;
252  }
253  ERES = stack->GetERes(z, a, EINC);
255 
256  SetEResAfterDetector(-1.);
257  return (EINC - ERES);
258  }
259 
260 protected:
261 
262  TString fFName;
263  KVList* fCalibrators;
264  KVList* fParticles;
265  KVList fAbsorbers;
266 
267  Double_t ELossActive(Double_t* x, Double_t* par);
268  Double_t EResDet(Double_t* x, Double_t* par);
269  Double_t RangeDet(Double_t* x, Double_t* par);
270 
271  TF1* fELossF;
272  TF1* fEResF;
273  TF1* fRangeF;
274 
275  Double_t fEResforEinc;
276 
277  Bool_t fSimMode;
278  Bool_t fPresent;
279  Bool_t fDetecting;
280 
281  Bool_t fSingleLayer;
282 
283  void AddDetectorSignal(KVDetectorSignal* ds)
284  {
288  ds->SetDetector(this);
289  fDetSignals.Add(ds);
290  }
291 
292 public:
293  KVDetector();
294  KVDetector(const Char_t* type, const Float_t thick = 0.0);
295  KVDetector(const Char_t* gas, const Double_t thick, const Double_t pressure, const Double_t temperature = 19.0);
296  KVDetector(const KVDetector&);
297  void init();
298  virtual ~ KVDetector();
299 
300  void Copy(TObject& obj) const override;
301 
302  void SetMaterial(const Char_t* type) override;
303  void AddAbsorber(KVMaterial*);
304  void SetActiveLayer(KVMaterial* actif)
305  {
306  fActiveLayer = actif;
307  SetBit(kActiveSet);
308  }
309  void SetActiveLayer(Int_t i)
310  {
312  if (auto mat = GetAbsorber(i)) SetActiveLayer(mat);
313  }
314  KVMaterial* GetActiveLayer() const override
315  {
317  return fActiveLayer;
318  }
319  KVMaterial* GetAbsorber(Int_t i) const;
320  KVMaterial* GetAbsorber(const Char_t* name) const
321  {
323  return (KVMaterial*)fAbsorbers.FindObject(name);
324  }
325  const KVList* GetListOfAbsorbers() const
326  {
327  return &fAbsorbers;
328  }
329  Int_t GetNumberOfAbsorberLayers() const
330  {
331  return fAbsorbers.GetEntries();
332  }
333  void RemoveAllAbsorbers();
334 
335  Double_t GetTotalThicknessInCM() const
336  {
339 
340  double fTotThickness = 0;
341  TIter next(&fAbsorbers);
342  KVMaterial* mat;
343  while ((mat = (KVMaterial*)next())) fTotThickness += mat->GetThickness();
344  return fTotThickness;
345  }
346  KVGeoDetectorNode* GetNode()
347  {
348  return &fNode;
349  }
350 
351  const Char_t* GetMaterialName() const
352  {
353  if (GetActiveLayer())
354  return GetActiveLayer()->GetName();
355  return KVMaterial::GetName();
356  }
357  void DetectParticle(KVNucleus*, TVector3* norm = 0) override;
358  Double_t GetELostByParticle(KVNucleus*, TVector3* norm = 0) override;
359  Double_t GetParticleEIncFromERes(KVNucleus*, TVector3* norm = 0) override;
360 
361  virtual Double_t GetCalibratedEnergy() const
362  {
364  return GetDetectorSignalValue("Energy");
365  }
366  virtual Double_t GetEnergy() const
367  {
371  if (IsSimMode()) return ELoss; // in simulation mode, return calculated energy loss in active layer
372  if (ELoss > 0) return ELoss;
373  ELoss = GetCalibratedEnergy();
374  if (ELoss < 0) ELoss = 0;
375  SetEnergy(ELoss);
376  return ELoss;
377  }
378  virtual void SetEnergy(Double_t e) const
379  {
385  }
386  Double_t GetEnergyLoss() const override
387  {
388  return GetEnergy();
389  }
390  void SetEnergyLoss(Double_t e) const override
391  {
392  SetEnergy(e);
393  }
394  virtual Double_t GetCorrectedEnergy(KVNucleus*, Double_t e =
395  -1., Bool_t transmission = kTRUE);
396  virtual Int_t FindZmin(Double_t ELOSS = -1., Char_t mass_formula = -1);
397 
398  Bool_t AddCalibrator(KVCalibrator* cal, const KVNameValueList& opts = "");
399  Bool_t ReplaceCalibrator(const Char_t* type, KVCalibrator* cal, const KVNameValueList& opts = "");
400  KVCalibrator* GetCalibrator(const Char_t* name,
401  const Char_t* type) const;
402  KVCalibrator* GetCalibrator(const Char_t* type) const;
403  KVList* GetListOfCalibrators() const
404  {
405  return fCalibrators;
406  }
407  Bool_t IsCalibrated() const
408  {
410  return (IsSimMode() && IsDetecting()) || (HasDetectorSignal("Energy") && IsOK());
411  }
412  Bool_t IsCalibrated(const KVNameValueList& params) const;
413 
414  void Clear(Option_t* opt = "") override;
415  virtual void Reset(Option_t* opt = "")
416  {
417  Clear(opt);
418  }
419  void Print(Option_t* option = "") const override;
420 
421  void AddHit(KVNucleus* part)
422  {
424 
425  if (!fParticles) {
426  fParticles = new KVList(kFALSE);
427  fParticles->SetCleanup();
428  }
429  fParticles->Add(part);
430  }
431 
432  void RemoveHit(KVNucleus* part)
433  {
435 
436  fParticles->Remove(part);
437  if (!fParticles->GetEntries()) SetAnalysed(kFALSE);
438  }
439 
441  KVList* GetHits() const
442  {
443  return fParticles;
444  }
445  void ClearHits()
446  {
448  if (fParticles) fParticles->Clear();
449  }
451  Int_t GetNHits() const
452  {
453  return (fParticles ? fParticles->GetEntries() : 0);
454  }
455 
456  Bool_t IsAnalysed()
457  {
458  return TestBit(kIsAnalysed);
459  }
460  void SetAnalysed(Bool_t b = kTRUE)
461  {
462  SetBit(kIsAnalysed, b);
463  }
464  virtual Bool_t Fired(Option_t* opt = "any") const
465  {
480 
481  if (!IsDetecting()) return kFALSE; //detector not working, no answer at all
482  if (IsSimMode()) return (GetActiveLayer() ? GetActiveLayer()->GetEnergyLoss() > 0. : KVMaterial::GetEnergyLoss() > 0.); // simulation mode: detector fired if energy lost in active layer
483 
484  TString OPT(opt);
485  OPT.ToLower();
486  Bool_t all = (OPT == "all");
487 
488  TIter raw_it(&GetListOfDetectorSignals());
489  KVDetectorSignal* ds;
490  int count_raw = 0;
491  while ((ds = (KVDetectorSignal*)raw_it())) {
492  if (ds->IsRaw()) {
493  ++count_raw;
494  if (ds->IsFired()) {
495  if (!all) return kTRUE;
496  }
497  else {
498  if (all) return kFALSE;
499  }
500  }
501  }
502  return all && count_raw;
503  }
504  virtual void RemoveCalibrators();
505 
506  Double_t GetDetectorSignalValue(const KVString& type, const KVNameValueList& params = "") const
507  {
515 
516  KVDetectorSignal* s = GetDetectorSignal(type);
517  return (s ? s->GetValue(params) : 0);
518  }
519  void SetDetectorSignalValue(const KVString& type, Double_t val) const
520  {
526 
527  KVDetectorSignal* s = GetDetectorSignal(type);
528  if (s) s->SetValue(val);
529  }
530  Double_t GetInverseDetectorSignalValue(const KVString& output, Double_t value, const KVString& input, const KVNameValueList& params = "") const
531  {
542 
543  KVDetectorSignal* s = GetDetectorSignal(output);
544  return (s ? s->GetInverseValue(value, input, params) : 0);
545  }
546  virtual KVDetectorSignal* GetDetectorSignal(const KVString& type) const
547  {
552 
553  return fDetSignals.get_object<KVDetectorSignal>(type);
554  }
555  Bool_t HasDetectorSignal(const KVString& type) const
556  {
559  return (GetDetectorSignal(type) != nullptr);
560  }
561 
562  inline void IncrementUnidentifiedParticles(Int_t n = 1)
563  {
564  fUnidentP += n;
565  fUnidentP = (fUnidentP > 0) * fUnidentP;
566  SetBit(kUnidentifiedParticle, (Bool_t)(fUnidentP > 0));
567  }
568  inline void IncrementIdentifiedParticles(Int_t n = 1)
569  {
570  fIdentP += n;
571  fIdentP = (fIdentP > 0) * fIdentP;
572  SetBit(kIdentifiedParticle, (Bool_t)(fIdentP > 0));
573  }
574  Bool_t BelongsToUnidentifiedParticle() const
575  {
576  return TestBit(kUnidentifiedParticle);
577  }
578  Bool_t BelongsToIdentifiedParticle() const
579  {
580  return TestBit(kIdentifiedParticle);
581  }
582 
583  static KVDetector* MakeDetector(const Char_t* name, Float_t thick);
584 
585  virtual Double_t GetEntranceWindowSurfaceArea();
586  TVector3 GetActiveLayerSurfaceCentre() const
587  {
590  return GetSurfaceCentre();
591  }
592  TVector3 GetActiveLayerVolumeCentre() const
593  {
596  }
597  TGeoBBox* GetActiveLayerShape() const
598  {
600  return GetShape();
601  }
602  TGeoHMatrix* GetActiveLayerMatrix() const
603  {
605  return GetMatrix();
606  }
607 
608  Double_t GetMaxDeltaE(Int_t Z, Int_t A) override;
609  Double_t GetEIncOfMaxDeltaE(Int_t Z, Int_t A) override;
610  Double_t GetDeltaE(Int_t Z, Int_t A, Double_t Einc, Double_t = 0.) override;
611  virtual Double_t GetTotalDeltaE(Int_t Z, Int_t A, Double_t Einc);
612  Double_t GetERes(Int_t Z, Int_t A, Double_t Einc, Double_t = 0.) override;
613  Double_t GetIncidentEnergy(Int_t Z, Int_t A, Double_t delta_e =
614  -1.0, enum SolType type = kEmax) override;
615  /*virtual Double_t GetEResFromDeltaE(...) - DON'T IMPLEMENT, CALLS GETINCIDENTENERGY*/
616  Double_t GetDeltaEFromERes(Int_t Z, Int_t A, Double_t Eres) override;
617  Double_t GetIncidentEnergyFromERes(Int_t Z, Int_t A, Double_t Eres) override;
618  Double_t GetRange(Int_t Z, Int_t A, Double_t Einc) override;
619  Double_t GetLinearRange(Int_t Z, Int_t A, Double_t Einc) override;
620  Double_t GetPunchThroughEnergy(Int_t Z, Int_t A) override;
621  virtual TGraph* DrawPunchThroughEnergyVsZ(Int_t massform = KVNucleus::kBetaMass);
622  virtual TGraph* DrawPunchThroughEsurAVsZ(Int_t massform = KVNucleus::kBetaMass);
623 
624  virtual TF1* GetEResFunction(Int_t Z, Int_t A);
625  virtual TF1* GetELossFunction(Int_t Z, Int_t A);
626  virtual TF1* GetRangeFunction(Int_t Z, Int_t A);
627 
628  virtual Double_t GetSmallestEmaxValid(Int_t Z, Int_t A) const;
629 
630  virtual void SetEResAfterDetector(Double_t e)
631  {
632  fEResforEinc = e;
633  }
634  virtual Double_t GetEResAfterDetector() const
635  {
636  return fEResforEinc;
637  }
638 
639  virtual void ReadDefinitionFromFile(const Char_t*);
640 
641  virtual void SetSimMode(Bool_t on = kTRUE)
642  {
648  fSimMode = on;
649  }
650  virtual Bool_t IsSimMode() const
651  {
657  return fSimMode;
658  }
659 
660  virtual Bool_t IsPresent() const
661  {
663  return fPresent;
664  }
665  void SetPresent(Bool_t yes = kTRUE)
666  {
667  fPresent = yes;
668  }
669  virtual Bool_t IsDetecting() const
670  {
672  return fDetecting;
673  }
674  void SetDetecting(Bool_t yes = kTRUE)
675  {
676  fDetecting = yes;
677  }
678 
679  virtual Bool_t IsOK() const
680  {
682  return (fPresent && fDetecting);
683  }
684 
685  KVGroup* GetGroup() const;
686  UInt_t GetGroupNumber();
687 
688  void AddParentStructure(KVGeoStrucElement* elem);
689  void RemoveParentStructure(KVGeoStrucElement* elem);
690  KVGeoStrucElement* GetParentStructure(const Char_t* type, const Char_t* name = "") const;
691 
692  void SetActiveLayerMatrix(const TGeoHMatrix*);
693  void SetActiveLayerShape(TGeoBBox*);
694  void SetEntranceWindowMatrix(const TGeoHMatrix*);
695  void SetEntranceWindowShape(TGeoBBox*);
696  const KVPosition& GetEntranceWindow() const
697  {
700  return fEWPosition;
701  }
702  const TVector3 GetCentreOfEntranceWindow() const
703  {
707  return GetEntranceWindow().GetSurfaceCentre();
708  }
709  Double_t GetSolidAngle() const override
710  {
712  if (ROOTGeo()) return fEWPosition.GetSolidAngle();
713  return KVPosition::GetSolidAngle();
714  }
715  TVector3 GetRandomDirection(Option_t* t = "isotropic") override
716  {
718  if (ROOTGeo()) return fEWPosition.GetRandomDirection(t);
720  }
721  void GetRandomAngles(Double_t& th, Double_t& ph, Option_t* t = "isotropic") override
722  {
724  if (ROOTGeo()) fEWPosition.GetRandomAngles(th, ph, t);
725  else KVPosition::GetRandomAngles(th, ph, t);
726  }
727  TVector3 GetDirection() override
728  {
730  if (ROOTGeo()) return fEWPosition.GetDirection();
731  return KVPosition::GetDirection();
732  }
733  Double_t GetDistance() const override
734  {
736  if (ROOTGeo()) return fEWPosition.GetDistance();
737  return KVPosition::GetDistance();
738  }
739  Double_t GetTheta() const override
740  {
742  if (ROOTGeo()) return fEWPosition.GetTheta();
743  return KVPosition::GetTheta();
744  }
745  Double_t GetSinTheta() const override
746  {
748  if (ROOTGeo()) return fEWPosition.GetSinTheta();
749  return KVPosition::GetSinTheta();
750  }
751  Double_t GetCosTheta() const override
752  {
754  if (ROOTGeo()) return fEWPosition.GetCosTheta();
755  return KVPosition::GetCosTheta();
756  }
757  Double_t GetPhi() const override
758  {
760  if (ROOTGeo()) return fEWPosition.GetPhi();
761  return KVPosition::GetPhi();
762  }
763 
764  void SetThickness(Double_t thick) override;
765  Bool_t IsSingleLayer() const
766  {
768  return fSingleLayer;
769  }
770  Bool_t IsMultiLayer() const
771  {
773  return !IsSingleLayer();
774  }
775  Bool_t IsGasDetector() const
776  {
778  return IsGas();
779  }
780  Bool_t HasSameStructureAs(const KVDetector*) const;
781  void SetNameOfArray(const TString& n)
782  {
783  fNameOfArray = n;
784  }
785  const Char_t* GetNameOfArray() const
786  {
788  return fNameOfArray;
789  }
790 
791  const KVSeqCollection& GetListOfDetectorSignals() const
792  {
793  return fDetSignals;
794  }
795  KVDetectorSignal* AddDetectorSignal(const KVString& type)
796  {
801 
802  auto signal = new KVDetectorSignal(type, this);
803  fDetSignals.Add(signal);
804  return signal;
805  }
806  Bool_t AddDetectorSignalExpression(const KVString& type, const KVString& _expr);
807 
808  virtual Int_t GetIndex() const
809  {
813 
814  return 0;
815  }
816 
817  void SetPressure(Double_t P) override;
818  void SetTemperature(Double_t T) override;
819 
820  virtual Bool_t IsSegmented() const
821  {
822  return kFALSE;
823  }
824  virtual void AddEnergyLossInSubDetector(int, double) {}
825 
826  ClassDefOverride(KVDetector, 10) //Base class for the description of detectors in multidetector arrays
827 };
828 
829 inline KVCalibrator* KVDetector::GetCalibrator(const Char_t* name,
830  const Char_t* type) const
831 {
832  if (fCalibrators)
833  return (KVCalibrator*) fCalibrators->FindObjectWithNameAndType(name, type);
834  return 0;
835 }
836 
837 inline KVCalibrator* KVDetector::GetCalibrator(const Char_t* type) const
838 {
839  if (fCalibrators)
840  return (KVCalibrator*) fCalibrators->FindObjectByType(type);
841  return 0;
842 }
843 #endif
int Int_t
unsigned int UInt_t
#define OPT
#define e(i)
bool Bool_t
char Char_t
float Float_t
constexpr Bool_t kFALSE
double Double_t
constexpr Bool_t kTRUE
const char Option_t
#define BIT(n)
#define ClassDefOverride(name, id)
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void on
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t Int_t Int_t Window_t TString Int_t GCValues_t GetPrimarySelectionOwner GetDisplay GetScreen GetColormap GetNativeEvent const char const char dpyName wid window const char font_name cursor keysym reg const char only_if_exist regb h Point_t winding char text const char depth char const char Int_t count const char ColorStruct_t color const char Pixmap_t Pixmap_t PictureAttributes_t attr const char char ret_data h unsigned char height h Atom_t Int_t ULong_t ULong_t unsigned char prop_list Atom_t Atom_t Atom_t Time_t type
Base class for all detector calibrations.
Definition: KVCalibrator.h:99
Base class for output signal data produced by a detector.
virtual Bool_t IsRaw() const
Bool_t IsFired() const
void SetDetector(const KVDetector *d)
Information on relative positions of detectors & particle trajectories.
Base class describing elements of array geometry.
Group of detectors which can be treated independently of all others in array.
Definition: KVGroup.h:20
Extended TList class which owns its objects by default.
Definition: KVList.h:28
Description of physical materials used to construct detectors & targets; interface to range tables.
Definition: KVMaterial.h:89
virtual void SetEnergyLoss(Double_t e) const
Definition: KVMaterial.h:149
virtual Double_t GetThickness() const
Definition: KVMaterial.cpp:484
virtual Double_t GetEnergyLoss() const
Definition: KVMaterial.h:139
virtual Double_t GetEIncOfMaxDeltaE(Int_t Z, Int_t A)
void Clear(Option_t *opt="") override
Reset absorber - set stored energy lost by particles in absorber to zero.
void init()
Definition: KVMaterial.cpp:46
virtual Double_t GetParticleEIncFromERes(KVNucleus *, TVector3 *norm=nullptr)
Definition: KVMaterial.cpp:801
virtual void DetectParticle(KVNucleus *, TVector3 *norm=nullptr)
virtual void SetMaterial(const Char_t *type)
Definition: KVMaterial.cpp:217
virtual Double_t GetELostByParticle(KVNucleus *, TVector3 *norm=nullptr)
Definition: KVMaterial.cpp:766
void Copy(TObject &obj) const override
Make a copy of this material object.
void Print(Option_t *option="") const override
Show information on this material.
Definition: KVMaterial.cpp:739
virtual KVMaterial * GetActiveLayer() const
Definition: KVMaterial.h:176
Handles lists of named parameters with different types, a list of KVNamedParameter objects.
Description of properties and kinematics of atomic nuclei.
Definition: KVNucleus.h:123
Int_t GetA() const
Definition: KVNucleus.cpp:792
Int_t GetZ() const
Return the number of proton / atomic number.
Definition: KVNucleus.cpp:763
Base class used for handling geometry in a multidetector array.
Definition: KVPosition.h:91
virtual void GetRandomAngles(Double_t &th, Double_t &ph, Option_t *t="isotropic")
Definition: KVPosition.cpp:280
virtual TVector3 GetRandomDirection(Option_t *t="isotropic")
Definition: KVPosition.cpp:242
virtual void SetShape(TGeoBBox *)
Definition: KVPosition.cpp:743
virtual Double_t GetSolidAngle(void) const
Definition: KVPosition.cpp:567
virtual Double_t GetTheta() const
Definition: KVPosition.h:160
virtual TVector3 GetSurfaceCentre() const
Definition: KVPosition.cpp:886
virtual TGeoHMatrix * GetMatrix() const
Definition: KVPosition.cpp:777
virtual Double_t GetPhi() const
Definition: KVPosition.h:172
virtual Double_t GetDistance(void) const
Definition: KVPosition.h:190
virtual TVector3 GetVolumeCentre() const
Definition: KVPosition.cpp:915
virtual Double_t GetSurfaceArea(int npoints=100000) const
Definition: KVPosition.cpp:974
virtual Double_t GetSinTheta() const
Definition: KVPosition.h:164
virtual TVector3 GetSurfaceNormal() const
Definition: KVPosition.cpp:945
virtual Double_t GetMisalignmentAngle() const
virtual void SetMatrix(const TGeoHMatrix *)
Definition: KVPosition.cpp:702
virtual TGeoBBox * GetShape() const
Definition: KVPosition.cpp:790
virtual Double_t GetCosTheta() const
Definition: KVPosition.h:168
virtual TVector3 GetDirection()
Definition: KVPosition.cpp:463
virtual TVector3 GetRandomPointOnSurface() const
Definition: KVPosition.cpp:813
KaliVeda extensions to ROOT collection classes.
T * get_object(const TString &name) const
TObject * Remove(TObject *obj) override
Remove object from list.
void Add(TObject *obj) override
TObject * FindObject(const char *name) const override
void Clear(Option_t *option="") override
virtual void SetCleanup(Bool_t enable=kTRUE)
virtual TObject * FindObjectByType(const Char_t *) const
virtual TObject * FindObjectWithNameAndType(const Char_t *name, const Char_t *type) const
Extension of ROOT TString class which allows backwards compatibility with ROOT v3....
Definition: KVString.h:73
Optimised list in which named objects can only be placed once.
void Add(TObject *obj) override
virtual Int_t GetEntries() const
const char * GetName() const override
void SetBit(UInt_t f)
R__ALWAYS_INLINE Bool_t TestBit(UInt_t f) const
const Int_t n
const long double s
Definition: KVUnits.h:94
TArc a