MDetector/Module.cc

Go to the documentation of this file.

#include <mdet/Module.h>

#include <mdet/Counter.h>
#include <mdet/MHierarchyInfo.h>
#include <sdet/Station.h>
#include <utl/Vector.h>
#include <utl/ErrorLogger.h>
#include <vector>
#include <sstream>

#include <boost/iterator/transform_iterator.hpp>
#include <boost/lambda/lambda.hpp>


namespace mdet {

  const char* const Module::kComponentName = MHierarchyInfo::kComponentsNames[1];


  const char* const Module::kComponentId = MHierarchyInfo::kComponentsIds[1];


  const char* const Module::AreaKindTags[] = { "small", "large" };

  const char* const Module::PhotoDetectorTags[] = { "PMT", "SiPM" };

  const Counter&
  Module::GetCounter()
    const
  {
    return fCounter;
  }


  const PMT&
  Module::GetPMT()
    const
  {
    return *fPMT;
  }

  const SiPMArray&
  Module::GetSiPMArray()
    const
  {
    return *fSiPMArray;
  }

  const FrontEnd&
  Module::GetFrontEnd()
    const
  {
    return *fFrontEnd;
  }

  const FrontEndSiPM&
  Module::GetFrontEndSiPM()
    const
  {
    return *fFrontEndSiPM;
  }

  const BackEndSiPM&
  Module::GetBackEndSiPM()
    const
  {
    return *fBackEndSiPM;
  }

  Module::Module(const int mId, const det::VManager::IndexMap& parentMap, const Counter& parent) :
    MDetectorComponent<Module>::Type(mId, parentMap),
    det::MPositionable<Module>(*this),
    fScintillators(*this),
    fFibers(*this),
    fCounter(parent),
    fIndexByChannel(fIndexRelations.get<ChannelTag>()),
    fIndexByPixel(fIndexRelations.get<PixelTag>()),
    fIndexByScintillator(fIndexRelations.get<ScintillatorTag>()),
    fIndexByFiber(fIndexRelations.get<FiberTag>())
  {
    fScintillators.Update(GetIdsMap());
    fFibers.Update(GetIdsMap());
    // Being only one, the id is fixed to one.
    unsigned int numberOfPixels = 0;
    unsigned int numberOfChannels = 0;
    if(GetPhotoDetector()==0){
        fPMT = new PMT(1, GetIdsMap(), *this);
        fFrontEnd = new FrontEnd(1, GetIdsMap(), *this);
        numberOfPixels = fPMT->GetPixels().GetNumberOfComponents();
        numberOfChannels = fFrontEnd->GetChannels().GetNumberOfComponents();
    }else{
        fSiPMArray = new SiPMArray(1, GetIdsMap(), *this);
        fFrontEndSiPM = new FrontEndSiPM(1, GetIdsMap(), *this);
        fBackEndSiPM = new BackEndSiPM(1, GetIdsMap(), *this);
        numberOfPixels = fSiPMArray->GetSiPMs().GetNumberOfComponents();
        numberOfChannels = fFrontEndSiPM->GetChannels().GetNumberOfComponents();
    }
    // Check length of vectors.
    if (!(fScintillators.GetNumberOfComponents() == fFibers.GetNumberOfComponents() &&
          fScintillators.GetNumberOfComponents() == numberOfPixels &&
          fScintillators.GetNumberOfComponents() == numberOfChannels)) {

      if(GetPhotoDetector()==0)
         INFO("This Module uses PMT.");
      else
         INFO("This Module uses SiPM.");

      std::ostringstream ss;
      ss << "For the module ";
      AddIdMessage(ss);
      ss << " there are "
         << fScintillators.GetNumberOfComponents() << " scintillators, "
         << fFibers.GetNumberOfComponents() << " fibers, "
         << numberOfPixels << " pixels/sipms, "
         << numberOfChannels << " channels "
            " when the quantities must be the same.";
      FATAL(ss);
      throw;
    }
    Register(fAreaKind);
    Register(fBackgroundMuonsFlux);
    Register(fPhotoDetector);
  }


  utl::CoordinateSystemPtr
  Module::GetReferenceCoordinateSystem()
    const
  {
    return GetCounter().GetLocalCoordinateSystem();
  }


  Module::AreaKind
  Module::GetAreaKind()
    const
  {
    return AreaKindCreator::Create(GetData(fAreaKind, "areaKind"));
  }


  double
  Module::GetBackgroundMuonsFlux()
    const
  {
    return GetData(fBackgroundMuonsFlux, "backgroundMuonsFlux");
  }


  Module::PhotoDetector
  Module::GetPhotoDetector()
    const
  {
    return PhotoDetectorCreator::Create(GetData(fPhotoDetector, "photoDetector"));
  }

  double
  Module::GetPhi0()
    const
  {
    return GetEulerPhi();
  }


  double
  Module::GetArea()
    const
  {
    double area = 0;

    for (ScintillatorConstIterator is = ScintillatorsBegin(); is != ScintillatorsEnd(); ++is)
      area += is->GetArea();

    return area;
  }


  unsigned int
  Module::GetNumberOfScintillators()
    const
  {
    unsigned int nScintillators = 0;

    for (ScintillatorConstIterator is = ScintillatorsBegin(); is != ScintillatorsEnd(); ++is)
      ++nScintillators;

    return nScintillators;
  }


  ULong64_t
  Module::GetMask()
    const
  {
    return fFrontEnd->GetMask();
  }


  double
  Module::GetDepth()
    const
  {
    return GetDepth(GetPosition());
  }


  double
  Module::GetDepth(const utl::Point& p)
    const
  {
    /*
     * What's better: to couple or not to the station/counter
     * this depth? At first seems to be the natural reference;
     * and it wouldn't make much sense to have an arbitrary depth.
     *
     * TODO Check if sdet::GetPosition effectively returns the on-
     * the-ground position of the tank, with X-Y defining the ground.
     */
    const sdet::Station& station = GetCounter().GetAssociatedTank();
    // Construct the distance vector, ending in the point, and starting
    // in the (surface) station position
    const utl::Vector dist = p - station.GetPosition();
    // Use -1 because we want this depth to be positive for an underground point.
    const utl::Vector normal2ground(0, 0, -1, station.GetLocalCoordinateSystem());
    // Project the Z component: should be positive for underground modules (will
    // there be any "flying" module...? hahaah)
    return dist * normal2ground;
  }

  bool
  Module::IsSiPM()
    const
  {
    return (bool) (GetPhotoDetector() == eSiPM);
  }

  double
  Module::GetShortestFiber()
    const
  {
          double length = 5000;
          for (Module::FiberConstIterator fIt = FibersBegin(); fIt != FibersEnd(); ++fIt) {
                  length = std::min(length, fIt->GetOnManifoldLength());
          }
          return length;
  }

  void
  Module::Update(const bool invalidateData, const bool invalidateComponents)
  {
    MDetectorComponent<Module>::Type::Update(invalidateData, invalidateComponents);
    fScintillators.Update(GetIdsMap(), invalidateData, invalidateComponents);
    fFibers.Update(GetIdsMap(), invalidateData, invalidateComponents);

    INFO("rereading channel-pmt/sipm-scintillator assignment..");
    fIndexByChannel.clear();
    if(!IsSiPM()){
                for (mdet::PMT::PixelConstIterator pixIt = fPMT->PixelsBegin(); pixIt != fPMT->PixelsEnd(); ++pixIt) {
                  const Pixel& pix = *pixIt;
                  const IdsContainer ids = pix.GetChannelScintillatorFiberIds();
                  if (ids.size() != 3) {
                        std::ostringstream ss;
                        ss << "The information in MModelConfig is corrupt! Number of indexes in <channelScintillatorFiberIds> has to be three";
                        FATAL(ss);
                        throw;
                  }


                  IdTuple newtuple;
                  newtuple.fScintillatorId = ids[1];
                  newtuple.fFiberId = ids[2];
                  newtuple.fPixelId = pixIt->GetId();
                  newtuple.fChannelId = ids[0];
                  fIndexByChannel.insert(newtuple);
                }
    }else{
                for (mdet::SiPMArray::SiPMConstIterator sipmIt = fSiPMArray->SiPMsBegin(); sipmIt != fSiPMArray->SiPMsEnd(); ++sipmIt) {
                  const SiPM& sipm = *sipmIt;
                  const IdsContainer ids = sipm.GetChannelSiPMScintillatorFiberIds();
                  if (ids.size() != 3) {
                        std::ostringstream ss;
                        ss << "The information in MModelConfig is corrupt! Number of indexes in <channelScintillatorFiberIds> has to be three";
                        FATAL(ss);
                        throw;
                  }

                  IdTuple newtuple;
                  newtuple.fScintillatorId = ids[1];
                  newtuple.fFiberId = ids[2];
                  newtuple.fPixelId = sipmIt->GetId();
                  newtuple.fChannelId = ids[0];
                  fIndexByChannel.insert(newtuple);
                }

    }
  }
}