KVPosition.h

/***************************************************************************
                          kvposition.h  -  description
                             -------------------
    begin                : Sun May 19 2002
    copyright            : (C) 2002 by J.D. Frankland
    email                : frankland@ganil.fr
 
$Id: KVPosition.h,v 1.20 2009/01/19 14:36:57 ebonnet Exp $
 ***************************************************************************/
 
/***************************************************************************
 *                                                                         *
 *   This program is free software; you can redistribute it and/or modify  *
 *   it under the terms of the GNU General Public License as published by  *
 *   the Free Software Foundation; either version 2 of the License, or     *
 *   (at your option) any later version.                                   *
 *                                                                         *
 ***************************************************************************/
 
#ifndef KVPOSITION_H
#define KVPOSITION_H
#include "TMath.h"
#include "KVBase.h"
#include "TVector3.h"
#include "TRotation.h"
#include "TGeoMatrix.h"
#include "TGeoBBox.h"
 
class KVPosition
{
private:
   Double_t fTheta = -1;              // polar angle in degrees with respect to beam axis, corresponds to centre of telescope
   mutable Double_t fPhi = -1;                // azimuthal angle in degrees with respect to 12 o'clock (=0 deg.), corresponds to centre of telescope
   Double_t fTheta_min = 0;          // polar angle in degrees corresponding to edge of telescope closest to beam axis
   Double_t fTheta_max= 0;          // polar angle in degrees of the edge furthest from the beam axis
   Double_t fPhi_min = 0;            // azimuthal angle in degrees corresponding to most anticlockwise edge of telescope
   Double_t fPhi_max = 0;            // azimuthal angle in degrees corresponding to most clockwise edge of telescope
   Double_t fDistance = 0;           // distance in cm from centre of solid angle element to coordinate system origin (target)
 
   /* ROOT Geometry */
    std::unique_ptr<TGeoHMatrix> fMatrix;
   TGeoBBox*    fShape = nullptr;         // shape of detector volume
   mutable Double_t     fSolidAngle = 0;  // solid angle = area of entrance window / distance**2
    TVector3 fVolumeCentre;
    TVector3 fSurfaceCentre;
    TVector3 fNormal;
   Double_t fMaxMisalignment=-1.;             // maximum angle [deg] between the normal & all vectors which pass through the surface
 
   Double_t GetOC_SC_CosAngle() const
   {
      return GetVolumeCentre().Unit().Dot(GetNormalIntoShape());
   }
   Double_t GetMisalignmentAngle() const
   {
       return TMath::RadToDeg() * TMath::ACos(GetOC_SC_CosAngle());
   }
 
   void initialize_root_geo();
   void set_volume_centre();
   void set_surface_centre();
 
   public:
   KVPosition() = default;
   KVPosition(Double_t thmin, Double_t thmax, Double_t phmin, Double_t phmax, Double_t dist = 0.0);
   virtual ~KVPosition() = default;
   KVPosition(const KVPosition& o)
       : fTheta(o.fTheta),fPhi(o.fPhi),fTheta_min(o.fTheta_min),
       fTheta_max(o.fTheta_max),fPhi_min(o.fPhi_min),fDistance(o.fDistance),
       fMatrix(o.GetMatrix() ? new TGeoHMatrix(*o.GetMatrix()) : nullptr),
       fShape(o.fShape),fSolidAngle(o.fSolidAngle),
       fVolumeCentre(o.fVolumeCentre),
       fSurfaceCentre(o.fSurfaceCentre),
       fNormal(o.fNormal) {}
   KVPosition(KVPosition&& other) noexcept
       : KVPosition()
   {
       swap(*this,other);
   }
 
   friend void swap(KVPosition& p1, KVPosition& p2)
   {
       using std::swap;
       swap(p1.fTheta,p2.fTheta);
       swap(p1.fPhi,p2.fPhi);
       swap(p1.fTheta_min,p2.fTheta_min);
       swap(p1.fTheta_max,p2.fTheta_max);
       swap(p1.fPhi_min,p2.fPhi_min);
       swap(p1.fPhi_max,p2.fPhi_max);
       swap(p1.fDistance,p2.fDistance);
       swap(p1.fMatrix,p2.fMatrix);
       swap(p1.fShape,p2.fShape);
       swap(p1.fSolidAngle,p2.fSolidAngle);
       swap(p1.fVolumeCentre,p2.fVolumeCentre);
       swap(p1.fSurfaceCentre,p2.fSurfaceCentre);
       swap(p1.fNormal,p2.fNormal);
   }
   KVPosition& operator=(KVPosition other)
   {
       swap(*this, other);
       return *this;
   }
 
   virtual void SetAzimuthalAngle(Double_t ph);
   virtual void SetPolarAngle(Double_t th);
   virtual void SetPolarWidth(Double_t pw);
   virtual void SetPolarMinMax(Double_t min, Double_t max);
   virtual void SetAzimuthalWidth(Double_t aw);
   virtual void SetAzimuthalMinMax(Double_t min, Double_t max);
   void SetPhiMinMax(Double_t min, Double_t max)
   {
      SetAzimuthalMinMax(min, max);
   }
   virtual TVector3 GetRandomDirection(Option_t* t = "isotropic");
   TRotation GetRandomIsotropicRotation();
   virtual void GetRandomAngles(Double_t& th, Double_t& ph, Option_t* t = "isotropic");
   void CalculateEstimatedAngularLimits(int N=10000);
   virtual TVector3 GetDirection();
 
   Bool_t IsInPolarRange(const Double_t theta);
   Bool_t IsInPhiRange(const Double_t phi);
   Bool_t IsSmallerThan(KVPosition* pos);
   Bool_t IsAlignedWith(KVPosition* pos);
   Bool_t IsOverlappingWith(KVPosition* pos);
   Bool_t IsAzimuthallyWiderThan(KVPosition* pos);
 
   Double_t GetPhiMin() const
   {
      if(ROOTGeo() && fPhi_min==fPhi_max)
         const_cast<KVPosition*>(this)->CalculateEstimatedAngularLimits();
      return fPhi_min;
   }
   Double_t GetPhiMax() const
   {
      if(ROOTGeo() && fPhi_min==fPhi_max)
         const_cast<KVPosition*>(this)->CalculateEstimatedAngularLimits();
      return fPhi_max;
   }
   Double_t GetThetaMin() const
   {
      if(ROOTGeo() && fTheta_min==fTheta_max)
         const_cast<KVPosition*>(this)->CalculateEstimatedAngularLimits();
      return fTheta_min;
   }
   Double_t GetThetaMax() const
   {
      if(ROOTGeo() && fTheta_min==fTheta_max)
         const_cast<KVPosition*>(this)->CalculateEstimatedAngularLimits();
      return fTheta_max;
   }
   Double_t GetAzimuthalWidth(Double_t phmin = -1., Double_t phimax =
                                 -1.) const;
   virtual Double_t GetTheta() const
   {
      return fTheta;
   }
   virtual Double_t GetSinTheta() const
   {
      return TMath::Sin(fTheta * TMath::DegToRad());
   }
   virtual Double_t GetCosTheta() const
   {
      return TMath::Cos(fTheta * TMath::DegToRad());
   }
   virtual Double_t GetPhi() const
   {
      if (fPhi >= 360.0)
         fPhi -= 360.;
      return fPhi;
   }
   void SetTheta(Double_t t)
   {
      SetPolarAngle(t);
   }
   void SetPhi(Double_t p)
   {
      SetAzimuthalAngle(p);
   }
   void SetDistance(Double_t d)
   {
      fDistance = d;
   }
   virtual Double_t GetDistance(void) const
   {
      return fDistance;
   }
   void GetCornerCoordinates(TVector3*, Double_t /*depth*/ = 0);
   void GetCornerCoordinatesInOwnFrame(TVector3*, Double_t /*depth*/ = 0);
   void GetWidthsFromDimension(Double_t lin_dim);
 
   virtual Double_t GetSolidAngle(void) const;
 
   /* ROOT Geometry */
   virtual Bool_t ROOTGeo() const
   {
      return (fMatrix && fShape);
   }
   virtual void SetMatrix(const TGeoHMatrix*);
   virtual void SetShape(TGeoBBox*);
   virtual TGeoHMatrix* GetMatrix() const;
   virtual TGeoBBox* GetShape() const;
   virtual TVector3 GetRandomPointOnSurface() const;
   virtual TVector3 GetSurfaceCentre() const
   {
       return fSurfaceCentre;
   }
   virtual TVector3 GetVolumeCentre() const
   {
       return fVolumeCentre;
   }
   virtual TVector3 GetNormalIntoShape() const
   {
       return fNormal;
   }
   virtual Double_t GetSurfaceArea(int npoints = 100000) const;
   virtual Double_t GetMaximumMisalignmentAngle() const
   {
 
       if(fMaxMisalignment<0)
           const_cast<KVPosition*>(this)->CalculateEstimatedAngularLimits();
       return fMaxMisalignment;
   }
   virtual Double_t GetCosineMaximumMisalignmentAngle() const
   {
 
       return TMath::Cos(GetMaximumMisalignmentAngle()*TMath::DegToRad());
   }
 
   ClassDef(KVPosition, 2)//Class handling geometry of detectors in a multidetector array
};
 
 
 
#endif