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  KVUniqueNameList fDetSignalsForRawTree;
223 
224  void remove_signal_for_calibrator(KVCalibrator* K);
225 
226  template<typename AbsorberStack>
227  Double_t get_corrected_energy(AbsorberStack* stack, KVNucleus* nuc, Double_t e, Bool_t transmission)
228  {
229  auto z = nuc->GetZ();
230  auto a = nuc->GetA();
231  enum SolType solution = kEmax;
232  if (!transmission) solution = kEmin;
233  Double_t EINC, ERES = GetEResAfterDetector();
234  if (transmission && ERES > 0.) { // if residual energy is known we use it to calculate EINC.
235 
236  EINC = stack->GetIncidentEnergyFromERes(z, a, ERES);
237  if (EINC < stack->GetEIncOfMaxDeltaE(z, a)) {// if EINC < max of dE curve, we change solution
238  solution = kEmin;
239  } // we could keep the EINC value calculated using ERES, but then the corrected dE of this detector would not depend on the measured dE !
240  }
241  EINC = stack->GetIncidentEnergy(z, a, e, solution);
242  if (EINC < 0) {
243  SetEResAfterDetector(-1.);
244  return EINC;
245  }
246  ERES = stack->GetERes(z, a, EINC);
247 
248  SetEResAfterDetector(-1.);
249  return (EINC - ERES);
250  }
251 
252 protected:
253 
254  TString fFName;
255  KVList* fCalibrators;
256  KVList* fParticles;
257  KVList fAbsorbers;
258 
259  Double_t ELossActive(Double_t* x, Double_t* par);
260  Double_t EResDet(Double_t* x, Double_t* par);
261  Double_t RangeDet(Double_t* x, Double_t* par);
262 
263  TF1* fELossF;
264  TF1* fEResF;
265  TF1* fRangeF;
266 
267  Double_t fEResforEinc;
268 
269  Bool_t fSimMode;
270  Bool_t fPresent;
271  Bool_t fDetecting;
272 
273  Bool_t fSingleLayer;
274 
275  void AddDetectorSignal(KVDetectorSignal* ds)
276  {
280  ds->SetDetector(this);
281  fDetSignals.Add(ds);
282  if (use_signal_for_raw_data_tree(ds->GetType())) fDetSignalsForRawTree.Add(ds);
283  }
284  virtual Bool_t use_signal_for_raw_data_tree(const TString&) const
285  {
286  return true;
287  }
288 
289 public:
290  KVDetector();
291  KVDetector(const Char_t* type, const Float_t thick = 0.0);
292  KVDetector(const Char_t* gas, const Double_t thick, const Double_t pressure, const Double_t temperature = 19.0);
293  KVDetector(const KVDetector&);
294  void init();
295  virtual ~ KVDetector();
296 
297  void Copy(TObject& obj) const override;
298 
299  void SetMaterial(const Char_t* type) override;
300  void AddAbsorber(KVMaterial*);
301  void SetActiveLayer(KVMaterial* actif)
302  {
303  fActiveLayer = actif;
304  SetBit(kActiveSet);
305  }
306  void SetActiveLayer(Int_t i)
307  {
309  if (auto mat = GetAbsorber(i)) SetActiveLayer(mat);
310  }
311  KVMaterial* GetActiveLayer() const override
312  {
314  return fActiveLayer;
315  }
316  KVMaterial* GetAbsorber(Int_t i) const;
317  KVMaterial* GetAbsorber(const Char_t* name) const
318  {
320  return (KVMaterial*)fAbsorbers.FindObject(name);
321  }
322  const KVList* GetListOfAbsorbers() const
323  {
324  return &fAbsorbers;
325  }
326  Int_t GetNumberOfAbsorberLayers() const
327  {
328  return fAbsorbers.GetEntries();
329  }
330  void RemoveAllAbsorbers();
331 
332  Double_t GetTotalThicknessInCM() const
333  {
336 
337  double fTotThickness = 0;
338  TIter next(&fAbsorbers);
339  KVMaterial* mat;
340  while ((mat = (KVMaterial*)next())) fTotThickness += mat->GetThickness();
341  return fTotThickness;
342  }
343  KVGeoDetectorNode* GetNode()
344  {
345  return &fNode;
346  }
347 
348  const Char_t* GetMaterialName() const
349  {
350  if (GetActiveLayer())
351  return GetActiveLayer()->GetName();
352  return KVMaterial::GetName();
353  }
354  void DetectParticle(KVNucleus*, TVector3* norm = 0) override;
355  Double_t GetELostByParticle(KVNucleus*, TVector3* norm = 0) override;
356  Double_t GetParticleEIncFromERes(KVNucleus*, TVector3* norm = 0) override;
357 
358  virtual Double_t GetCalibratedEnergy() const
359  {
361  return GetDetectorSignalValue("Energy");
362  }
363  virtual Double_t GetEnergy() const
364  {
368  if (IsSimMode()) return ELoss; // in simulation mode, return calculated energy loss in active layer
369  if (ELoss > 0) return ELoss;
370  ELoss = GetCalibratedEnergy();
371  if (ELoss < 0) ELoss = 0;
372  SetEnergy(ELoss);
373  return ELoss;
374  }
375  virtual void SetEnergy(Double_t e) const
376  {
382  }
383  Double_t GetEnergyLoss() const override
384  {
385  return GetEnergy();
386  }
387  void SetEnergyLoss(Double_t e) const override
388  {
389  SetEnergy(e);
390  }
391  virtual Double_t GetCorrectedEnergy(KVNucleus*, Double_t e =
392  -1., Bool_t transmission = kTRUE);
393  virtual Int_t FindZmin(Double_t ELOSS = -1., Char_t mass_formula = -1);
394 
395  Bool_t AddCalibrator(KVCalibrator* cal, const KVNameValueList& opts = "");
396  Bool_t ReplaceCalibrator(const Char_t* type, KVCalibrator* cal, const KVNameValueList& opts = "");
397  KVCalibrator* GetCalibrator(const Char_t* name,
398  const Char_t* type) const;
399  KVCalibrator* GetCalibrator(const Char_t* type) const;
400  KVList* GetListOfCalibrators() const
401  {
402  return fCalibrators;
403  }
404  Bool_t IsCalibrated() const
405  {
407  return (IsSimMode() && IsDetecting()) || (HasDetectorSignal("Energy") && IsOK());
408  }
409  Bool_t IsCalibrated(const KVNameValueList& params) const;
410 
411  void Clear(Option_t* opt = "") override;
412  virtual void Reset(Option_t* opt = "")
413  {
414  Clear(opt);
415  }
416  void Print(Option_t* option = "") const override;
417 
418  void AddHit(KVNucleus* part)
419  {
421 
422  if (!fParticles) {
423  fParticles = new KVList(kFALSE);
424  fParticles->SetCleanup();
425  }
426  fParticles->Add(part);
427  }
428 
429  void RemoveHit(KVNucleus* part)
430  {
432 
433  fParticles->Remove(part);
434  if (!fParticles->GetEntries()) SetAnalysed(kFALSE);
435  }
436 
438  KVList* GetHits() const
439  {
440  return fParticles;
441  }
442  void ClearHits()
443  {
445  if (fParticles) fParticles->Clear();
446  }
448  Int_t GetNHits() const
449  {
450  return (fParticles ? fParticles->GetEntries() : 0);
451  }
452 
453  Bool_t IsAnalysed()
454  {
455  return TestBit(kIsAnalysed);
456  }
457  void SetAnalysed(Bool_t b = kTRUE)
458  {
459  SetBit(kIsAnalysed, b);
460  }
461  virtual Bool_t Fired(Option_t* opt = "any") const
462  {
477 
478  if (!IsDetecting()) return kFALSE; //detector not working, no answer at all
479  if (IsSimMode()) return (GetActiveLayer() ? GetActiveLayer()->GetEnergyLoss() > 0. : KVMaterial::GetEnergyLoss() > 0.); // simulation mode: detector fired if energy lost in active layer
480 
481  TString OPT(opt);
482  OPT.ToLower();
483  Bool_t all = (OPT == "all");
484 
485  TIter raw_it(&GetListOfDetectorSignals());
486  KVDetectorSignal* ds;
487  int count_raw = 0;
488  while ((ds = (KVDetectorSignal*)raw_it())) {
489  if (ds->IsRaw()) {
490  ++count_raw;
491  if (ds->IsFired()) {
492  if (!all) return kTRUE;
493  }
494  else {
495  if (all) return kFALSE;
496  }
497  }
498  }
499  return all && count_raw;
500  }
501  virtual void RemoveCalibrators();
502 
503  Double_t GetDetectorSignalValue(const KVString& type, const KVNameValueList& params = "") const
504  {
512 
513  KVDetectorSignal* s = GetDetectorSignal(type);
514  return (s ? s->GetValue(params) : 0);
515  }
516  void SetDetectorSignalValue(const KVString& type, Double_t val) const
517  {
523 
524  KVDetectorSignal* s = GetDetectorSignal(type);
525  if (s) s->SetValue(val);
526  }
527  Double_t GetInverseDetectorSignalValue(const KVString& output, Double_t value, const KVString& input, const KVNameValueList& params = "") const
528  {
539 
540  KVDetectorSignal* s = GetDetectorSignal(output);
541  return (s ? s->GetInverseValue(value, input, params) : 0);
542  }
543  virtual KVDetectorSignal* GetDetectorSignal(const KVString& type) const
544  {
549 
550  return fDetSignals.get_object<KVDetectorSignal>(type);
551  }
552  Bool_t HasDetectorSignal(const KVString& type) const
553  {
556  return (GetDetectorSignal(type) != nullptr);
557  }
558 
559  inline void IncrementUnidentifiedParticles(Int_t n = 1)
560  {
561  fUnidentP += n;
562  fUnidentP = (fUnidentP > 0) * fUnidentP;
563  SetBit(kUnidentifiedParticle, (Bool_t)(fUnidentP > 0));
564  }
565  inline void IncrementIdentifiedParticles(Int_t n = 1)
566  {
567  fIdentP += n;
568  fIdentP = (fIdentP > 0) * fIdentP;
569  SetBit(kIdentifiedParticle, (Bool_t)(fIdentP > 0));
570  }
571  Bool_t BelongsToUnidentifiedParticle() const
572  {
573  return TestBit(kUnidentifiedParticle);
574  }
575  Bool_t BelongsToIdentifiedParticle() const
576  {
577  return TestBit(kIdentifiedParticle);
578  }
579 
580  static KVDetector* MakeDetector(const Char_t* name, Float_t thick);
581 
582  virtual Double_t GetEntranceWindowSurfaceArea();
583  TVector3 GetActiveLayerSurfaceCentre() const
584  {
587  return GetSurfaceCentre();
588  }
589  TVector3 GetActiveLayerVolumeCentre() const
590  {
593  }
594  TGeoBBox* GetActiveLayerShape() const
595  {
597  return GetShape();
598  }
599  TGeoHMatrix* GetActiveLayerMatrix() const
600  {
602  return GetMatrix();
603  }
604 
605  Double_t GetMaxDeltaE(Int_t Z, Int_t A) override;
606  Double_t GetEIncOfMaxDeltaE(Int_t Z, Int_t A) override;
607  Double_t GetDeltaE(Int_t Z, Int_t A, Double_t Einc, Double_t = 0.) override;
608  virtual Double_t GetTotalDeltaE(Int_t Z, Int_t A, Double_t Einc);
609  Double_t GetERes(Int_t Z, Int_t A, Double_t Einc, Double_t = 0.) override;
610  Double_t GetIncidentEnergy(Int_t Z, Int_t A, Double_t delta_e =
611  -1.0, enum SolType type = kEmax) override;
612  /*virtual Double_t GetEResFromDeltaE(...) - DON'T IMPLEMENT, CALLS GETINCIDENTENERGY*/
613  Double_t GetDeltaEFromERes(Int_t Z, Int_t A, Double_t Eres) override;
614  Double_t GetIncidentEnergyFromERes(Int_t Z, Int_t A, Double_t Eres) override;
615  Double_t GetRange(Int_t Z, Int_t A, Double_t Einc) override;
616  Double_t GetLinearRange(Int_t Z, Int_t A, Double_t Einc) override;
617  Double_t GetPunchThroughEnergy(Int_t Z, Int_t A) override;
618  virtual TGraph* DrawPunchThroughEnergyVsZ(Int_t massform = KVNucleus::kBetaMass);
619  virtual TGraph* DrawPunchThroughEsurAVsZ(Int_t massform = KVNucleus::kBetaMass);
620 
621  virtual TF1* GetEResFunction(Int_t Z, Int_t A);
622  virtual TF1* GetELossFunction(Int_t Z, Int_t A);
623  virtual TF1* GetRangeFunction(Int_t Z, Int_t A);
624 
625  virtual Double_t GetSmallestEmaxValid(Int_t Z, Int_t A) const;
626 
627  virtual void SetEResAfterDetector(Double_t e)
628  {
629  fEResforEinc = e;
630  }
631  virtual Double_t GetEResAfterDetector() const
632  {
633  return fEResforEinc;
634  }
635 
636  virtual void ReadDefinitionFromFile(const Char_t*);
637 
638  virtual void SetSimMode(Bool_t on = kTRUE)
639  {
645  fSimMode = on;
646  }
647  virtual Bool_t IsSimMode() const
648  {
654  return fSimMode;
655  }
656 
657  virtual Bool_t IsPresent() const
658  {
660  return fPresent;
661  }
662  void SetPresent(Bool_t yes = kTRUE)
663  {
664  fPresent = yes;
665  }
666  virtual Bool_t IsDetecting() const
667  {
669  return fDetecting;
670  }
671  void SetDetecting(Bool_t yes = kTRUE)
672  {
673  fDetecting = yes;
674  }
675 
676  virtual Bool_t IsOK() const
677  {
679  return (fPresent && fDetecting);
680  }
681 
682  KVGroup* GetGroup() const;
683  UInt_t GetGroupNumber();
684 
685  void AddParentStructure(KVGeoStrucElement* elem);
686  void RemoveParentStructure(KVGeoStrucElement* elem);
687  KVGeoStrucElement* GetParentStructure(const Char_t* type, const Char_t* name = "") const;
688 
689  void SetActiveLayerMatrix(const TGeoHMatrix*);
690  void SetActiveLayerShape(TGeoBBox*);
691  void SetEntranceWindowMatrix(const TGeoHMatrix*);
692  void SetEntranceWindowShape(TGeoBBox*);
693  const KVPosition& GetEntranceWindow() const
694  {
697  return fEWPosition;
698  }
699  const TVector3 GetCentreOfEntranceWindow() const
700  {
704  return GetEntranceWindow().GetSurfaceCentre();
705  }
706  Double_t GetSolidAngle() const override
707  {
709  if (ROOTGeo()) return fEWPosition.GetSolidAngle();
710  return KVPosition::GetSolidAngle();
711  }
712  TVector3 GetRandomDirection(Option_t* t = "isotropic") override
713  {
715  if (ROOTGeo()) return fEWPosition.GetRandomDirection(t);
717  }
718  void GetRandomAngles(Double_t& th, Double_t& ph, Option_t* t = "isotropic") override
719  {
721  if (ROOTGeo()) fEWPosition.GetRandomAngles(th, ph, t);
722  else KVPosition::GetRandomAngles(th, ph, t);
723  }
724  TVector3 GetDirection() override
725  {
727  if (ROOTGeo()) return fEWPosition.GetDirection();
728  return KVPosition::GetDirection();
729  }
730  Double_t GetDistance() const override
731  {
733  if (ROOTGeo()) return fEWPosition.GetDistance();
734  return KVPosition::GetDistance();
735  }
736  Double_t GetTheta() const override
737  {
739  if (ROOTGeo()) return fEWPosition.GetTheta();
740  return KVPosition::GetTheta();
741  }
742  Double_t GetSinTheta() const override
743  {
745  if (ROOTGeo()) return fEWPosition.GetSinTheta();
746  return KVPosition::GetSinTheta();
747  }
748  Double_t GetCosTheta() const override
749  {
751  if (ROOTGeo()) return fEWPosition.GetCosTheta();
752  return KVPosition::GetCosTheta();
753  }
754  Double_t GetPhi() const override
755  {
757  if (ROOTGeo()) return fEWPosition.GetPhi();
758  return KVPosition::GetPhi();
759  }
760 
761  void SetThickness(Double_t thick) override;
762  Bool_t IsSingleLayer() const
763  {
765  return fSingleLayer;
766  }
767  Bool_t IsMultiLayer() const
768  {
770  return !IsSingleLayer();
771  }
772  Bool_t IsGasDetector() const
773  {
775  return IsGas();
776  }
777  Bool_t HasSameStructureAs(const KVDetector*) const;
778  void SetNameOfArray(const TString& n)
779  {
780  fNameOfArray = n;
781  }
782  const Char_t* GetNameOfArray() const
783  {
785  return fNameOfArray;
786  }
787 
788  const KVSeqCollection& GetListOfDetectorSignals() const
789  {
790  return fDetSignals;
791  }
792  const KVSeqCollection& GetListOfDetectorSignalsForRawTree() const
793  {
794  return fDetSignalsForRawTree;
795  }
796  KVDetectorSignal* AddDetectorSignal(const KVString& type)
797  {
802 
803  auto signal = new KVDetectorSignal(type, this);
804  fDetSignals.Add(signal);
805  if (use_signal_for_raw_data_tree(type)) fDetSignalsForRawTree.Add(signal);
806  return signal;
807  }
808  Bool_t AddDetectorSignalExpression(const KVString& type, const KVString& _expr);
809 
810  virtual Int_t GetIndex() const
811  {
815 
816  return 0;
817  }
818 
819  void SetPressure(Double_t P) override;
820  void SetTemperature(Double_t T) override;
821 
822  virtual Bool_t IsSegmented() const
823  {
824  return kFALSE;
825  }
826  virtual void AddEnergyLossInSubDetector(int, double) {}
827 
828  ClassDefOverride(KVDetector, 10) //Base class for the description of detectors in multidetector arrays
829 };
830 
831 inline KVCalibrator* KVDetector::GetCalibrator(const Char_t* name,
832  const Char_t* type) const
833 {
834  if (fCalibrators)
835  return (KVCalibrator*) fCalibrators->FindObjectWithNameAndType(name, type);
836  return 0;
837 }
838 
839 inline KVCalibrator* KVDetector::GetCalibrator(const Char_t* type) const
840 {
841  if (fCalibrators)
842  return (KVCalibrator*) fCalibrators->FindObjectByType(type);
843  return 0;
844 }
845 #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.
const Char_t * GetType() const override
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