RPCCosmicSeedrecHitFinder.cc

Go to the documentation of this file.

#include "RecoMuon/MuonSeedGenerator/src/RPCCosmicSeedrecHitFinder.h"
#include <DataFormats/TrackingRecHit/interface/TrackingRecHit.h>
#include <DataFormats/RPCRecHit/interface/RPCRecHit.h>
#include <FWCore/Framework/interface/ESHandle.h>
#include <Geometry/Records/interface/MuonGeometryRecord.h>

using namespace std;
using namespace edm;

MuonTransientTrackingRecHit::MuonRecHitContainer::const_iterator find(
    MuonTransientTrackingRecHit::MuonRecHitContainer::const_iterator firstIter,
    MuonTransientTrackingRecHit::MuonRecHitContainer::const_iterator lastIter,
    const MuonTransientTrackingRecHit::MuonRecHitPointer& recHitRef) {
  MuonTransientTrackingRecHit::MuonRecHitContainer::const_iterator index = lastIter;
  for (MuonTransientTrackingRecHit::MuonRecHitContainer::const_iterator iter = firstIter; iter != lastIter; iter++)
    if ((*iter) == recHitRef)
      index = iter;
  return index;
}

RPCCosmicSeedrecHitFinder::RPCCosmicSeedrecHitFinder() {
  // Initiate the member
  isLayerset = false;
  isConfigured = false;
  isInputset = false;
  isOutputset = false;
  isEdgeset = false;
  BxRange = 0;
  MaxDeltaPhi = 0;
  ClusterSet.clear();
  innerBounds.clear();
  isLayersmixed = false;
  LayersinRPC.clear();
  therecHits.clear();
  isOuterLayerfilled = false;
}

RPCCosmicSeedrecHitFinder::~RPCCosmicSeedrecHitFinder() {}

void RPCCosmicSeedrecHitFinder::configure(const edm::ParameterSet& iConfig) {
  // Set the configuration
  BxRange = iConfig.getParameter<unsigned int>("BxRange");
  MaxDeltaPhi = iConfig.getParameter<double>("MaxDeltaPhi");
  ClusterSet = iConfig.getParameter<std::vector<int> >("ClusterSet");

  // Set the signal open
  isConfigured = true;
}

void RPCCosmicSeedrecHitFinder::setInput(MuonRecHitContainer (&recHits)[RPCLayerNumber]) {
  for (unsigned int i = 0; i < RPCLayerNumber; i++) {
    AllrecHits[i].clear();
    AllrecHits[i] = recHits[i];
  }

  // Set the signal open
  isInputset = true;
}

void RPCCosmicSeedrecHitFinder::setEdge(const RPCGeometry& rpcGeometry) {
  // Find all chamber in RB1in and collect their surface
  const std::vector<DetId>& AllRPCId = rpcGeometry.detIds();
  for (std::vector<DetId>::const_iterator it = AllRPCId.begin(); it != AllRPCId.end(); it++) {
    RPCDetId RPCId(it->rawId());
    int Region = RPCId.region();
    int Station = RPCId.station();
    int Layer = RPCId.layer();
    if (Region == 0 && Station == 1 && Layer == 1) {
      const BoundPlane RPCChamberSurface = rpcGeometry.chamber(RPCId)->surface();
      innerBounds.push_back(RPCChamberSurface);
    }
  }

  // Set the signal open
  isEdgeset = true;
}

void RPCCosmicSeedrecHitFinder::unsetEdge() {
  // Clear all surfaces of chambers in RB1in
  innerBounds.clear();
  isEdgeset = false;
}

void RPCCosmicSeedrecHitFinder::unsetInput() {
  for (unsigned int i = 0; i < RPCLayerNumber; i++)
    AllrecHits[i].clear();
  isInputset = false;
}

void RPCCosmicSeedrecHitFinder::setOutput(RPCSeedFinder* Seed) {
  theSeed = Seed;
  // Set the signal open
  isOutputset = true;
}

void RPCCosmicSeedrecHitFinder::setLayers(const std::vector<unsigned int>& Layers) {
  LayersinRPC = Layers;
  // Set the signal open
  isLayerset = true;
}

void RPCCosmicSeedrecHitFinder::fillrecHits() {
  if (isLayerset == false || isConfigured == false || isOutputset == false || isInputset == false ||
      isEdgeset == false) {
    cout << "Not set the IO or not configured yet" << endl;
    return;
  }

  therecHits.clear();

  if (LayersinRPC.empty()) {
    cout << "Not set with any layers" << endl;
    LayersinRPC.clear();
    therecHits.clear();
    isLayerset = false;
  }

  // check the layers, 1=all barrel, 2=all endcap, 3=mix
  unsigned int Component = LayerComponent();
  if (Component == 3)
    isLayersmixed = true;
  else
    isLayersmixed = false;

  GlobalVector initVector(0, 0, 0);
  const MuonRecHitPointer recHitRef;
  isOuterLayerfilled = false;
  complete(initVector, recHitRef);

  // Unset the signal
  LayersinRPC.clear();
  therecHits.clear();
  isLayerset = false;
}

int RPCCosmicSeedrecHitFinder::LayerComponent() {
  bool isBarrel = false;
  bool isEndcap = false;
  for (std::vector<unsigned int>::const_iterator it = LayersinRPC.begin(); it != LayersinRPC.end(); it++) {
    if ((*it) < BarrelLayerNumber)
      isBarrel = true;
    if ((*it) >= BarrelLayerNumber && (*it) < (BarrelLayerNumber + EachEndcapLayerNumber * 2))
      isEndcap = true;
  }
  if (isBarrel == true && isEndcap == true)
    return 3;
  if (isEndcap == true)
    return 2;
  if (isBarrel == true)
    return 1;
  return 0;
}

bool RPCCosmicSeedrecHitFinder::complete(const GlobalVector& lastSegment, const MuonRecHitPointer& lastrecHitRef) {
  bool isrecHitsfound = false;

  for (unsigned int i = 0; i < RPCLayerNumber; i++)
    for (MuonRecHitContainer::const_iterator it = AllrecHits[i].begin(); it != AllrecHits[i].end(); it++) {
      cout << "Finding recHits from " << i << " th layer" << endl;
      // information for recHits
      GlobalPoint currentPosition = (*it)->globalPosition();
      int currentBX;

      // Check validation
      if (!(*it)->isValid())
        continue;

      // Check BX range, be sure there is only RPCRecHit in the MuonRecHitContainer when use the dynamic_cast
      TrackingRecHit* thisTrackingRecHit = (*it)->hit()->clone();
      // Should also delete the RPCRecHit object cast by dynamic_cast<> ?
      RPCRecHit* thisRPCRecHit = dynamic_cast<RPCRecHit*>(thisTrackingRecHit);
      currentBX = thisRPCRecHit->BunchX();
      int ClusterSize = thisRPCRecHit->clusterSize();
      delete thisTrackingRecHit;
      // Check BX
      if ((unsigned int)abs(currentBX) > BxRange)
        continue;

      // Check cluster size
      bool Clustercheck = false;
      if (ClusterSet.empty())
        Clustercheck = true;
      for (std::vector<int>::const_iterator CluIter = ClusterSet.begin(); CluIter != ClusterSet.end(); CluIter++)
        if (ClusterSize == (*CluIter))
          Clustercheck = true;
      if (Clustercheck != true)
        continue;

      cout << "Candidate recHit's position: " << currentPosition.x() << ", " << currentPosition.y() << ", "
           << currentPosition.z() << ". BX : " << currentBX << endl;
      // Fill 1st recHit from outer layers and rest recHits from all layers
      if (isOuterLayerfilled == false) {
        // Pick out the recHit from outer layers and fill it
        if (!isouterLayer(*it))
          continue;

        // If pass all, add to the seed
        GlobalVector currentSegment = GlobalVector(0, 0, 0);
        cout << "1st recHit's global position: " << currentPosition.x() << ", " << currentPosition.y() << ", "
             << currentPosition.z() << ". BX: " << currentBX << endl;
        isrecHitsfound = true;
        therecHits.push_back(*it);
        isOuterLayerfilled = true;
        complete(currentSegment, *it);
        // Remember to pop the recHit before add another one from the same layer!
        therecHits.pop_back();
        isOuterLayerfilled = false;
      } else {
        GlobalPoint lastPosition = lastrecHitRef->globalPosition();
        TrackingRecHit* lastTrackingRecHit = lastrecHitRef->hit()->clone();
        // Should also delete the RPCRecHit object cast by dynamic_cast<> ?
        RPCRecHit* lastRPCRecHit = dynamic_cast<RPCRecHit*>(lastTrackingRecHit);
        int lastBX = lastRPCRecHit->BunchX();
        delete lastTrackingRecHit;

        // Check the Y coordinate, shoule be lower than current one
        if (currentPosition.y() >= lastPosition.y())
          continue;

        // Check the BX, should be larger than current one
        if (currentBX < lastBX)
          continue;

        // If be the 2nd recHit, just fill it
        bool isinsideRegion = isinsideAngleRange(lastSegment, lastPosition, currentPosition);
        cout << "Check isinsideRegion: " << isinsideRegion << endl;
        if (!isinsideRegion)
          continue;

        // If cross the edge the recHit should belong to another seed
        bool iscrossanyEdge = iscorssEdge(lastrecHitRef, *it);
        cout << "Check iscrossanyEdge: " << iscrossanyEdge << endl;
        if (iscrossanyEdge)
          continue;

        // If pass all, add to the seed
        unsigned int NumberinSeed = therecHits.size();
        GlobalVector currentSegment = (GlobalVector)(currentPosition - lastPosition);
        cout << (NumberinSeed + 1) << "th recHit's global position: " << currentPosition.x() << ", "
             << currentPosition.y() << ", " << currentPosition.z() << ". BX: " << currentBX << endl;
        isrecHitsfound = true;
        therecHits.push_back(*it);

        // if could not find next recHit in the search path, and have enough recHits already, that is the candidate
        bool findNext = complete(currentSegment, *it);
        if (findNext == false && therecHits.size() > 3) {
          for (ConstMuonRecHitContainer::const_iterator iter = therecHits.begin(); iter != therecHits.end(); iter++)
            cout << "Find recHit in seed candidate : " << (*iter)->globalPosition().x() << ", "
                 << (*iter)->globalPosition().y() << ", " << (*iter)->globalPosition().z() << endl;
          checkandfill();
        }

        // Remember to pop the recHit before add another one from the same layer!
        therecHits.pop_back();
      }
    }
  return isrecHitsfound;
}

bool RPCCosmicSeedrecHitFinder::isouterLayer(const MuonRecHitPointer& recHitRef) {
  bool isinsideLayers = false;
  for (std::vector<unsigned int>::const_iterator it = LayersinRPC.begin(); it != LayersinRPC.end(); it++) {
    MuonRecHitContainer::const_iterator index = find(AllrecHits[*it].begin(), AllrecHits[*it].end(), recHitRef);
    if (index != AllrecHits[*it].end())
      isinsideLayers = true;
  }
  return isinsideLayers;
}

bool RPCCosmicSeedrecHitFinder::isinsideAngleRange(const GlobalVector& lastSegment,
                                                   const GlobalPoint& lastPosition,
                                                   const GlobalPoint& currentPosition) {
  bool isinsideAngle = true;
  GlobalVector SegVec = currentPosition - lastPosition;
  if (lastSegment.mag() != 0)
    if (fabs((lastSegment.phi() - SegVec.phi()).value()) > MaxDeltaPhi)
      isinsideAngle = false;

  return isinsideAngle;
}

bool RPCCosmicSeedrecHitFinder::iscorssEdge(const MuonRecHitPointer& lastrecHitRef,
                                            const MuonRecHitPointer& currentrecHitRef) {
  bool iscorss = false;

  // Check if 2 recHits corss the inner bounds
  GlobalPoint lastPosition = lastrecHitRef->globalPosition();
  GlobalPoint currentPosition = currentrecHitRef->globalPosition();
  GlobalPoint testPosition((lastPosition.x() + currentPosition.x()) / 2,
                           (lastPosition.y() + currentPosition.y()) / 2,
                           (lastPosition.z() + currentPosition.z()) / 2);
  /*
    for(std::vector<BoundPlane>::const_iterator it = innerBounds.begin(); it != innerBounds.end(); it++) {
        //SurfaceOrientation::Side TestSide0 = it->side(currentPosition, 0);
        //SurfaceOrientation::Side TestSide1 = it->side(lastPosition, 0);
        SurfaceOrientation::Side TestSide = it->side(testPosition, 0);
        //cout << "Side of currentPosition: " << TestSide0 << ", Side of lastPosition: " << TestSide1 << ", Side of middlePosition: " << TestSide << endl;
        //if(TestSide != SurfaceOrientation::positiveSide)
            //iscorss = true;
    }
    */

  // Check when mixLayer is not set
  if (isLayersmixed == false) {
    DetId lastId = lastrecHitRef->geographicalId();
    RPCDetId lastRPCId(lastId.rawId());
    int lastRegion = lastRPCId.region();
    DetId currentId = currentrecHitRef->geographicalId();
    RPCDetId currentRPCId(currentId.rawId());
    int currentRegion = currentRPCId.region();
    // Check if 2 recHits from different regions
    if (lastRegion != currentRegion)
      iscorss = true;
  }

  return iscorss;
}

void RPCCosmicSeedrecHitFinder::checkandfill() {
  if (therecHits.size() >= 3) {
    theSeed->setrecHits(therecHits);
    theSeed->seed();
  } else
    cout << "Layer less than 3, could not fill a RPCSeedFinder" << endl;
}