CellGeo.cxx

Go to the documentation of this file.

//////////////////////////////////////////////////////////////////////////
// \file    CellGeo.cxx
// \brief   Cell geometry construction and call routines
// \version $Id: CellGeo.cxx,v 1.7 2012-11-14 21:09:52 jcoelho Exp $
// \author  $Author: jcoelho $
//////////////////////////////////////////////////////////////////////////
#include <iostream>

#include "TGeoManager.h"
#include "TGeoBBox.h"
#include "TGeoMatrix.h"
#include "TGeoNode.h"
#include "TVector3.h"

#include "GeometryObjects/CellGeo.h"
#include "GeometryObjects/CellUniqueId.h"

namespace geo
{
  std::vector<TGeoHMatrix> CellGeo::fgParentMatrices;

  /// Construct a cell geometry
  CellGeo::CellGeo(const std::vector<const TGeoNode*>& n, int depth)
    : fTranslatedMatrix(0)
  {
    fCellNode = n[depth];

    // The current split between parentMatrix and the final step is a
    // compromise between speed and memory usage. This way, we can find our
    // transformation matrix from the world frame with only one multiply, but
    // we use something like 16x less storage than if we kept the final matrix
    // for every single cell.

    TGeoHMatrix parentMatrix = *n.front()->GetMatrix();
    // This loop is omitting the last entry in 'n'. That'll be available
    // through fCellNode.
    for(int i = 1; i < depth; ++i){
      parentMatrix.Multiply(n[i]->GetMatrix());
    }

    // We know that the geometry construction traverses depth-first, so if this
    // parent matrix has seen before it will actually have been immediately
    // prior.
    if(fgParentMatrices.empty() ||
       !(parentMatrix == fgParentMatrices.back())){
      fgParentMatrices.push_back(parentMatrix);

      // Every module in the detector is made of two halves, so count
      // 2*#modules to get the expected number of parent matrices.
      if(fgParentMatrices.size() == 64*14*12*2+1){
        std::cerr << "Warning: CellGeo::fgParentMatrices size exceeds "
                  << 64*14*12*2 << " expected." << std::endl;
      }

      // If all is well, this is the last matrix, so don't waste any additional
      // space in the vector.
      if(fgParentMatrices.size() == 64*14*12*2){
        fgParentMatrices.shrink_to_fit();
      }
    }
    // Either it matches the previous matrix, or we just added it to the list,
    // either way it's the last index.
    fParentMatrixIdx = fgParentMatrices.size()-1;

    // Construct the unique id for this cell
    fId = NodesToUniqueId(n, depth);
  }

  //......................................................................
  CellGeo::~CellGeo()
  {
    delete fTranslatedMatrix;
  }

  //......................................................................

  /// Transform a position from local frame to world frame
  /// \param local : 3D array. Position in the local frame  Input.
  /// \param world : 3D array. Position in the world frame. Returned.
  void CellGeo::LocalToWorld(const double* local, double* world) const
  {
    GetGeoMatrix().LocalToMaster(local,world);
  }

  //......................................................................

  /// Transform a 3-vector from local frame to world frame
  /// \param local : 3D array. Position in the local frame  Input.
  /// \param world : 3D array. Position in the world frame. Returned.
  void CellGeo::LocalToWorldVect(const double* local, double* world) const
  {
    GetGeoMatrix().LocalToMasterVect(local,world);
  }

  //......................................................................

  /// Transform a position from world frame to local frame
  /// \param world : 3D array. Position in the world frame. Input.
  /// \param local : 3D array. Position in the local frame  Returned.
  void CellGeo::WorldToLocal(const double* local, double* world) const
  {
    GetGeoMatrix().MasterToLocal(local,world);
  }

  //......................................................................

  /// Transform a 3-vector from world frame to local frame
  /// \param world : 3D array. Position in the world frame. Input.
  /// \param local : 3D array. Position in the local frame  Returned.
  void CellGeo::WorldToLocalVect(const double* local, double* world) const
  {
    GetGeoMatrix().MasterToLocalVect(local,world);
  }

  //......................................................................

  void CellGeo::TranslateCell(double dx, double dy, double dz)
  {

    /// Set rotation angles to 0 for now
    double phi = 0;
    double theta = 0;
    double psi = 0;

    /// Build the translation and rotation matrices
    /// phi, theta and psi are euler angles
    TGeoRotation *rot = new TGeoRotation("alignmentRot",phi,theta,psi);
    /// dx, dy and dz are the translations we want
    TGeoTranslation *trans = new TGeoTranslation("",dx,dy,dz);
    TGeoCombiTrans *comb = new TGeoCombiTrans(*trans,*rot);

    this->TranslateCell(comb);

  }

  //......................................................................

  void CellGeo::TranslateCell(TGeoCombiTrans* comb)
  {
    if(fTranslatedMatrix) delete fTranslatedMatrix;

    fTranslatedMatrix = new TGeoHMatrix(GetGeoMatrix());
    fTranslatedMatrix->MultiplyLeft(comb);
  }
  
  //...................................................................... 
   
  void CellGeo::RestoreOriginal()
  {
    delete fTranslatedMatrix;
  }

  //......................................................................

  /// Return the center position of a cell.
  /// \param xyz : 3-D array. The returned location.
  /// \param localz : Distance along the length of the cell
  /// (cm). Default is center of cell
  void CellGeo::GetCenter(double* xyz, double localz) const
  {
    static double xyzLocal[3] = {0,0,0};
    xyzLocal[2] = localz;
    this->LocalToWorld(xyzLocal, xyz);
  }

  /// Return the center position of a cell.
  /// \param xyz :a Tvector3 object. The returned location.
  /// \param localz : Distance along the length of the cell
  /// (cm). Default is center of cell
  void CellGeo::GetCenter(TVector3& xyz, double localz) const
  {
    double xyzArray[3];
    this->GetCenter(xyzArray, localz);
    xyz.SetXYZ(xyzArray[0], xyzArray[1], xyzArray[2]);
  }


  //......................................................................

  /// Return the distance from a point in the cell to the read out side
  /// of the cell
  /// \param localz : Z position in cm in local cell frame
  /// \returns Distance in cm
  double CellGeo::DistToReadOut(double localz) const
  {
    return this->HalfL()-localz;
  }

  //......................................................................

  double CellGeo::HalfW() const
  {
    return ((TGeoBBox*)fCellNode->GetVolume()->GetShape())->GetDX();
  }

  //......................................................................

  double CellGeo::HalfL() const
  {
    return ((TGeoBBox*)fCellNode->GetVolume()->GetShape())->GetDZ();
  }

  //......................................................................

  double CellGeo::HalfD() const
  {
    return ((TGeoBBox*)fCellNode->GetVolume()->GetShape())->GetDY();
  }

  //......................................................................
  TGeoHMatrix CellGeo::GetGeoMatrix() const
  {
    if(fTranslatedMatrix) return *fTranslatedMatrix;

    TGeoHMatrix ret = fgParentMatrices[fParentMatrixIdx];
    ret.Multiply(fCellNode->GetMatrix());
    return ret;
  }

} // end namespace geo
////////////////////////////////////////////////////////////////////////