CDetector/Station.cc

Go to the documentation of this file.

#include <iostream>
#include <sstream>
#include <boost/lexical_cast.hpp>

#include <det/Detector.h>
#include <cdet/Station.h>
#include <cdet/CDetector.h>
#include <cevt/Pad.h>
#include <sdet/Station.h>
#include <utl/ErrorLogger.h>
#include <utl/Math.h>
#include <utl/Point.h>
#include <utl/Line.h>
#include <utl/GeometryUtilities.h>
#include <utl/UTMPoint.h>
#include <utl/TimeStamp.h>
#include <utl/UTCDateTime.h>
#include <utl/TimeRange.h>
#include <utl/TimeInterval.h>
#include <utl/Reader.h>
#include <utl/AugerUnits.h>
#include <utl/CoordinateSystem.h>
#include <utl/AugerCoordinateSystem.h>
#include <utl/ReferenceEllipsoid.h>
#include <fwk/LocalCoordinateSystem.h>

using namespace det;
using namespace cdet;
using namespace std;
using namespace utl;

const double Station::kTimeBinSize = 1.0*utl::nanosecond;


Station::Station(const int id) :
  fId(id)
  ,fStationIdString(boost::lexical_cast<string>(id))
  //fPadVector
  //fPartnerIds
  //fAxes
  //fName
  //fInGrid
  //fRadius
  //fHeight
  //fThickness
  //fCrowns
  //fReferenceSystem
  //fPosition
  //fWaterAbsLength
  //fWaterRefractionIndex
  //fLinerReflectivity
  //fLinerSigmaAlpha
  //fLinerSpecularLobe
  //fLinerSpecularSpike
  //fGroupId
  //fCommissionRange
  ,fDense(false)
{
  fRPCPositions.clear();
  fRPCRotations.clear();
  fRPCCoordinateSystems.clear();
}


Station::~Station()
{
  fRPCPositions.clear();
  fRPCRotations.clear();
  fRPCCoordinateSystems.clear();
}


const string&
Station::GetName()
  const
{ 
  if (!fName.IsValid()) {
    GetStationData(fName.Get(), "name", "stationList", "station name");
    fName.SetValid();
  }
  return fName.Get();
}


const TimeRange&
Station::GetCommissionTimeRange()
  const
{
  if (!fCommissionTimeRange.IsValid()) {
    if (fDense) {  // Dense stations are managed by SDetector
      fCommissionTimeRange = 
        Detector::GetInstance().GetSDetector().GetStation(fId).GetCommissionTimeRange();
    }
    else {
      string timeString;
      GetStationData(timeString, "commission", "stationList", "commission time");
      const TimeStamp start =
        boost::lexical_cast<UTCDateTime>(timeString).GetTimeStamp();
      GetStationData(timeString, "decommission", "stationList", "commission time");
      const TimeStamp stop =
        boost::lexical_cast<UTCDateTime>(timeString).GetTimeStamp();
      fCommissionTimeRange = TimeRange(start, stop);
    }
  }
  return fCommissionTimeRange.Get();
}


Point
Station::GetPosition()
  const
{

  if (!fPosition.IsValid()) {
    if (fDense) {  // Dense stations are managed by SDetector
      fPosition = 
        Detector::GetInstance().GetSDetector().GetStation(fId).GetPosition();
    }
    else {
      double northing;
      GetStationData(northing, "northing", "stationList", "northing");
      double easting;
      GetStationData(easting, "easting", "stationList", "easting");
      double altitude;
      GetStationData(altitude, "altitude", "stationList", "altitude");
      
      string ellipsoid;
      GetStationData(ellipsoid, "ellipsoid", "stationList", "ellipsoid for UTM conversion");
      fEllipsoid = ReferenceEllipsoid::GetEllipsoidIDFromString(ellipsoid);
      string band;
      GetStationData(band, "band", "stationList", "band for UTM conversion");
      int zone;
      GetStationData(zone, "zone", "stationList", "zone for UTM conversion");
      
      const UTMPoint utmPosition(northing, easting, altitude, zone, band[0], fEllipsoid);
      const CoordinateSystemPtr cs = Detector::GetInstance().GetSiteCoordinateSystem();
      fPosition = utmPosition.GetPoint(cs);
    }
  }

  return fPosition.Get();
}


int
Station::GetNumberRPCChambers() 
  const
{
  if (!fNumberRPCChambers.IsValid()) {
    GetStationData(fNumberRPCChambers.Get(), "numberRPCChambers", "rpc_positions", "number RPC chambers");
    fNumberRPCChambers.SetValid();
  }
  return fNumberRPCChambers.Get();
}



const Point&
Station::GetRPCPosition(const unsigned int id)
  const
{
  if (fRPCPositions.find(id) == fRPCPositions.end()) {

    vector<double> coords;
    GetRPCData(boost::lexical_cast<string>(id),coords,"RPCPosition", "rpc_positions", "RPC coordinates");

    // Return the RPC chamber position in the SD station coordinate system, whose origin is at 
    // the center of the tank.

    const sdet::Station & station = Detector::GetInstance().GetSDetector().GetStation(fId);
    
    const CoordinateSystemPtr CS = station.GetLocalCoordinateSystem();

    const double halfHeight = station.GetHeight()/2.;
    const double thickness = station.GetThickness();
    
    fRPCPositions.insert(make_pair(id,Point(coords[0], coords[1], coords[2]-halfHeight-thickness, CS)));

  }

  return (fRPCPositions.find(id))->second;
}

double
Station::GetRPCRotation(const unsigned int id)
  const
{
  if (fRPCRotations.find(id) == fRPCRotations.end()) {

    double rotation = 0.0;
    GetRPCData(boost::lexical_cast<string>(id),rotation,"RPCRotation", "rpc_positions", "RPC rotation");

    // Return the RPC chamber rotation
    
    fRPCRotations.insert(make_pair(id,rotation)); 
    
  }

  return fRPCRotations.find(id)->second;
}

CoordinateSystemPtr
Station::GetRPCCoordinateSystem(const unsigned int id)
  const
{
  if (fRPCCoordinateSystems.find(id) == fRPCCoordinateSystems.end()) {
    const CoordinateSystemPtr & stationCS = 
      det::Detector::GetInstance().GetSDetector().GetStation(fId).GetLocalCoordinateSystem();

    const Point rpcCenter = GetRPCPosition(id) + 
      Vector(0.0,0.0,GetRPCHousingThickness()[2]+GetRPCHousingInnerDimensions()[2]/2.0,stationCS);
    
    const CoordinateSystemPtr CS = fwk::LocalCoordinateSystem::Create(rpcCenter);
    
    const double rot = GetRPCRotation(id);

    fRPCCoordinateSystems.insert(make_pair(id,CoordinateSystem::RotationZ(-rot, CS)));
  }
  
  return fRPCCoordinateSystems.find(id)->second;
}


vector<double>
Station::GetTankSupportTopSlabDimensions()
  const
{
  if (!fTankSupportTopSlabDimensions.IsValid()) {
    GetStationData(fTankSupportTopSlabDimensions.Get(), "top_slab", "tank_support", "top slab dimensions");
    fTankSupportTopSlabDimensions.SetValid();
  }
  return fTankSupportTopSlabDimensions.Get();
}


vector<double>
Station::GetTankSupportCentralFootDimensions()
  const
{
  if (!fTankSupportCentralFootDimensions.IsValid()) {
    GetStationData(fTankSupportCentralFootDimensions.Get(), "central_foot", "tank_support", "central foot dimensions");
    fTankSupportCentralFootDimensions.SetValid();
  }
  return fTankSupportCentralFootDimensions.Get();
}


vector<double>
Station::GetTankSupportCentralFootBaseDimensions()
  const
{
  if (!fTankSupportCentralFootBaseDimensions.IsValid()) {
    GetStationData(fTankSupportCentralFootBaseDimensions.Get(), "central_foot_base", "tank_support", "central foot base dimensions");
    fTankSupportCentralFootBaseDimensions.SetValid();
  }
  return fTankSupportCentralFootBaseDimensions.Get();
}


vector<double>
Station::GetTankSupportOuterFootDimensions()
  const
{
  if (!fTankSupportOuterFootDimensions.IsValid()) {
    GetStationData(fTankSupportOuterFootDimensions.Get(), "outer_foot", "tank_support", "outer foot dimensions");
    fTankSupportOuterFootDimensions.SetValid();
  }
  return fTankSupportOuterFootDimensions.Get();
}


vector<double>
Station::GetTankSupportOuterFootBaseDimensions()
  const
{
  if (!fTankSupportOuterFootBaseDimensions.IsValid()) {
    GetStationData(fTankSupportOuterFootBaseDimensions.Get(), "outer_foot_base", "tank_support", "outer foot base dimensions");
    fTankSupportOuterFootBaseDimensions.SetValid();
  }
  return fTankSupportOuterFootBaseDimensions.Get();
}


double
Station::GetTankSupportOuterFootDistanceToCenter()
  const
{
  if (!fTankSupportOuterFootDistanceToCenter.IsValid()) {
    GetStationData(fTankSupportOuterFootDistanceToCenter.Get(), "outer_foot_distance_to_center", "tank_support", "outer foot distance to center");
    fTankSupportOuterFootDistanceToCenter.SetValid();
  }
  return fTankSupportOuterFootDistanceToCenter.Get();
}


vector<double>
Station::GetRPCHousingInnerDimensions()
  const
{
  if (!fRPCHousingInnerDimensions.IsValid()) {
    GetStationData(fRPCHousingInnerDimensions.Get(), "inner_dimensions", "rpc_housing", "RPC housing inner dimensions");
    fRPCHousingInnerDimensions.SetValid();
  }
  return fRPCHousingInnerDimensions.Get();
}


vector<double>
Station::GetRPCHousingThickness()
  const
{
  if (!fRPCHousingThickness.IsValid()) {
    GetStationData(fRPCHousingThickness.Get(), "thickness", "rpc_housing", "RPC housing thickness");
    fRPCHousingThickness.SetValid();
  }
  return fRPCHousingThickness.Get();
}


double
Station::GetRPCSizeX()
  const
{
  if (!fRPCSizeX.IsValid()) {
    GetStationData(fRPCSizeX.Get(), "size_x", "rpc_chamber", "RPC size in X");
    fRPCSizeX.SetValid();
  }
  return fRPCSizeX.Get();
}


double
Station::GetRPCSizeY()
  const
{
  if (!fRPCSizeY.IsValid()) {
    GetStationData(fRPCSizeY.Get(), "size_y", "rpc_chamber", "RPC size in Y");
    fRPCSizeY.SetValid();
  }
  return fRPCSizeY.Get();
}


vector<unsigned int>
Station::GetListPadIds()
  const
{
  if (!fListPadIds.IsValid() || fListPadIds.Get().begin() == fListPadIds.Get().end()) {
    fListPadIds.SetValid();
    for (int ir=1; ir<GetNumberRPCChambers()+1;++ir)
      for (int row=1; row<GetNumberPadsY()+1 ;++row)
        for (int col=1; col<GetNumberPadsX()+1;++col)
          fListPadIds.Get().push_back(10000*ir+(row-1)*GetNumberPadsX() + col);
  }
  
  return fListPadIds.Get();
}


int 
Station::GetNumberPadsX() 
  const
{
  if (!fNumberPadsX.IsValid()) {
    GetStationData(fNumberPadsX.Get(), "number_pads_x", "rpc_readout", "Number of pads in X");
    fNumberPadsX.SetValid();
  }
  return fNumberPadsX.Get();
}    


int 
Station::GetNumberPadsY() 
  const
{
  if (!fNumberPadsY.IsValid()) {
    GetStationData(fNumberPadsY.Get(), "number_pads_y", "rpc_readout", "Number of pads in Y");
    fNumberPadsY.SetValid();
  }
  return fNumberPadsY.Get();
}    


double 
Station::GetPadSizeX() 
  const
{
  if (!fPadSizeX.IsValid()) {
    fPadSizeX.Get() = GetRPCSizeX()/float(GetNumberPadsX());
    fPadSizeX.SetValid();
  }
  return fPadSizeX.Get();
}


double 
Station::GetPadSizeY() 
  const
{
  if (!fPadSizeY.IsValid()) {
    fPadSizeY.Get() = GetRPCSizeY()/float(GetNumberPadsY());
    fPadSizeY.SetValid();
  }
  return fPadSizeY.Get();
}


double 
Station::GetPadArea() 
  const
{
  if (!fPadArea.IsValid()) {
    fPadArea.Get() = GetPadSizeX()*GetPadSizeY();
    fPadArea.SetValid();
  }
  return fPadArea.Get();
}


Point  
Station::GetPadPosition(const cevt::Pad & pad) 
  const
{
  return GetPadPosition(pad.GetId());
}


Point  
Station::GetPadPosition(const unsigned int & id) 
  const
{
  const unsigned int chamberId = int(id/10000);
  const unsigned int padId = id - chamberId*10000;  
    
  const double sizeX = GetRPCSizeX();
  const double sizeY = GetRPCSizeY();
  const unsigned int nPadsX = GetNumberPadsX();
  const double padSizeX = GetPadSizeX();
  const double padSizeY = GetPadSizeY();

  const unsigned int row = int((padId-1)/nPadsX)  + 1; 
  const unsigned int col = padId - int((row-1)*nPadsX);

  // Coords of cener of pad w.r.t. center of RPC chamber
  const float padX = -sizeX/2.0 + (col-1)*padSizeX + padSizeX/2.0; 
  const float padY = -sizeY/2.0 + (row-1)*padSizeY + padSizeY/2.0; 
  const float padZ = 0.0*meter; 

  // Define RPC Coordinate System centered in center of RPC
  const utl::CoordinateSystemPtr rpcCS = GetRPCCoordinateSystem(chamberId);

  return Point(padX,padY,padZ,rpcCS);
}


unsigned int   
Station::GetPadId(const Point & pos) 
  const
{
  const unsigned int nRPC = GetNumberRPCChambers();
  const unsigned int numPadsX = GetNumberPadsX();

  const float rpcSizeX = GetRPCSizeX();
  const float rpcSizeY = GetRPCSizeY();
  const float rpcSizeZ = GetRPCHousingInnerDimensions()[2];

  const float padSizeX = GetPadSizeX();
  const float padSizeY = GetPadSizeY();

  double posX = numeric_limits<double>::max();
  double posY = numeric_limits<double>::max();
  double posZ = numeric_limits<double>::max();

  unsigned int rpcId = 0;

  for (unsigned int i = 1; i < nRPC+1; ++i) {

    const CoordinateSystemPtr & CS = GetRPCCoordinateSystem(i);

    posX = pos.GetX(CS);
    posY = pos.GetY(CS);
    posZ = pos.GetZ(CS);

    if (abs(posX) < rpcSizeX/2.0 && abs(posY) < rpcSizeY/2.0 && abs(posZ) < rpcSizeZ/2.0) {
      rpcId = i;
      break;
    }
  }

  int padId = 0;
  
  if (abs(posX) < rpcSizeX/2.0 && abs(posY) < rpcSizeY/2.0) {
    int col = int((posX-(-rpcSizeX/2.0))/padSizeX) + 1;
    int row = int((posY-(-rpcSizeY/2.0))/padSizeY) + 1;
    
    if (row > 0 && col > 0)
      padId = (row-1)*numPadsX + col;
  }

  if (rpcId > 0 && padId > 0)
    padId = 10000*rpcId + padId;
  
  return padId;
}


double
Station::GetDeadTime() 
  const
{
  if (!fDeadTime.IsValid()) {
    GetStationData(fDeadTime.Get(), "dead_time", "rpc_readout", "Dead time");
    fDeadTime.SetValid();
  }
  return fDeadTime.Get();
}  


bool
Station::IsInGrid(const CDetectorConstants::GridIndex index)
  const
{
  if (!fInGrid.IsValid()) {
    InGridType& inGrid = fInGrid.Get();
    inGrid.reset();
    vector<bool> data;

    if (fDense)
      data.clear();
    else 
      GetStationData(data, "inGrid", "stationList", "get is-in-grid");
    
    for (size_t i = 0; i < data.size() && i < inGrid.size(); ++i)
      inGrid[i] = data[i];
    fInGrid.SetValid();
  }

  const InGridType& inGrid = fInGrid.Get();
  if (index == CDetectorConstants::eAny)
    return inGrid.any();
  else
    return inGrid[index];
}


void
Station::NotFoundAndThrow(const string& msg)
  const
{
  ostringstream err;
  err << "Did not find requested component: '" << msg
      << "' for station = " << fStationIdString;
  throw NonExistentComponentException(err.str());
}


void
Station::Update()
  const
{
  /*
   In this method we clean all the variables that might change with time and deleted whatever was allocated in the get*'s
   Data members that are NOT reset are:
     fName
     fInGrid
     fRadius
     fHeight
     fThickness
     fCommissionTime
     fDecommissionTime
     fPosition (unless it is a dense station)
     fStationIdString
     fPartnerIdsSet
     fPartnerIds
     fAxes
   */

  if (fDense) {
    fPosition.SetValid(false);
    fReferenceSystem = CoordinateSystemPtr(); // clear the ptr so CS can be reassigned
    fRPCPositions.clear();
    fRPCRotations.clear();
    fRPCCoordinateSystems.clear();
    fListPadIds.Get().clear();
  }
}

// Configure (x)emacs for this file ...
// Local Variables:
// mode: c++
// End: