EcalTBReadout.cc

Go to the documentation of this file.

#include <algorithm>
#include <iostream>

#include "FWCore/MessageLogger/interface/MessageLogger.h"
#include "FWCore/Utilities/interface/Exception.h"
#include "SimCalorimetry/EcalTestBeamAlgos/interface/EcalTBReadout.h"
#include "SimDataFormats/EcalTestBeam/interface/PEcalTBInfo.h"

EcalTBReadout::EcalTBReadout(const std::string theEcalTBInfoLabel) : ecalTBInfoLabel_(theEcalTBInfoLabel) {
  theTargetCrystal_ = -1;
  theTTlist_.reserve(1);
}

void EcalTBReadout::findTTlist(const int &crysId, const EcalTrigTowerConstituentsMap &etmap) {
  // search for the TT involved in the NCRYMATRIX x NCRYMATRIX
  // around the target crystal if a new target, otherwise use
  // the list already filled

  if (crysId == theTargetCrystal_) {
    return;
  }

  theTTlist_.clear();


  EBDetId theTargetId;
  std::vector<DetId>::const_iterator idItr = theDetIds->begin();
  unsigned int ncount = 0;
  bool found = false;

  while ((ncount < theDetIds->size()) && !found) {
    EBDetId thisEBdetid(idItr->rawId());
    if (thisEBdetid.ic() == crysId) {
      theTargetId = thisEBdetid;
      found = true;
    }
    ++idItr;
    ++ncount;
  }
  if (!found) {
    throw cms::Exception("ObjectNotFound", "Ecal TB target crystal not found in geometry");
    return;
  }
  theTargetCrystal_ = theTargetId.ic();


  int myEta = theTargetId.ieta();
  int myPhi = theTargetId.iphi();

  for (int icrysEta = (myEta - (NCRYMATRIX - 1) / 2); icrysEta <= (myEta + (NCRYMATRIX - 1) / 2); ++icrysEta) {
    for (int icrysPhi = (myPhi - (NCRYMATRIX - 1) / 2); icrysPhi <= (myPhi + (NCRYMATRIX - 1) / 2); ++icrysPhi) {

      EBDetId thisEBdetid;

      idItr = theDetIds->begin();
      ncount = 0;
      found = false;

      while ((ncount < theDetIds->size()) && !found) {
        EBDetId myEBdetid(idItr->rawId());
        if ((myEBdetid.ieta() == icrysEta) && (myEBdetid.iphi() == icrysPhi)) {
          thisEBdetid = myEBdetid;
          found = true;
        }
        ++idItr;
        ++ncount;
      }

      if (found) {
        EcalTrigTowerDetId thisTTdetId = etmap.towerOf(thisEBdetid);

        LogDebug("EcalDigi") << "Crystal to be readout: sequential id = " << thisEBdetid.ic() << " eta = " << icrysEta
                             << " phi = " << icrysPhi << " from TT = " << thisTTdetId;

        if (theTTlist_.empty() || (theTTlist_.size() == 1 && theTTlist_[0] != thisTTdetId)) {
          theTTlist_.push_back(thisTTdetId);
        } else {
          std::vector<EcalTrigTowerDetId>::iterator ttFound = find(theTTlist_.begin(), theTTlist_.end(), thisTTdetId);
          if (theTTlist_.size() > 1 && ttFound == theTTlist_.end() && *(theTTlist_.end()) != thisTTdetId) {
            theTTlist_.push_back(thisTTdetId);
          }
        }
      }
    }
  }

  edm::LogInfo("EcalDigi") << " TT to be read: ";
  for (unsigned int i = 0; i < theTTlist_.size(); ++i) {
    edm::LogInfo("EcalDigi") << " TT " << i << " " << theTTlist_[i];
  }
}

void EcalTBReadout::readOut(const EBDigiCollection &input,
                            EBDigiCollection &output,
                            const EcalTrigTowerConstituentsMap &etmap) {
  /*
  for(EBDigiCollection::const_iterator digiItr = input.begin();
      digiItr != input.end(); ++digiItr)
    {
      EcalTrigTowerDetId thisTTdetId=etmap.towerOf(digiItr->id());
      std::vector<EcalTrigTowerDetId>::iterator ttFound =
  find(theTTlist_.begin(), theTTlist_.end(), thisTTdetId); if ((ttFound !=
  theTTlist_.end()) || *(theTTlist_.end()) == thisTTdetId) {
        output.push_back(*digiItr);
      }
    }
  edm::LogInfo("EcalDigi") << "Read EB Digis: " << output.size();
  */

  std::cout << "%%%%%%% In readOut(), size=" << input.size() << std::endl;

  for (unsigned int digis = 0; digis < input.size(); ++digis) {
    EBDataFrame ebdf = input[digis];

    EcalTrigTowerDetId thisTTdetId = etmap.towerOf(ebdf.id());
    std::vector<EcalTrigTowerDetId>::iterator ttFound = find(theTTlist_.begin(), theTTlist_.end(), thisTTdetId);

    if ((ttFound != theTTlist_.end()) || *(theTTlist_.end()) == thisTTdetId) {
      output.push_back(ebdf.id());
      EBDataFrame ebdf2(output.back());
      std::copy(ebdf.frame().begin(), ebdf.frame().end(), ebdf2.frame().begin());
    }
  }
}

void EcalTBReadout::readOut(const EEDigiCollection &input,
                            EEDigiCollection &output,
                            const EcalTrigTowerConstituentsMap &etmap) {
  for (unsigned int digis = 0; digis < input.size(); ++digis) {
    EEDataFrame eedf(input[digis]);

    EcalTrigTowerDetId thisTTdetId(etmap.towerOf(eedf.id()));

    std::vector<EcalTrigTowerDetId>::iterator ttFound(find(theTTlist_.begin(), theTTlist_.end(), thisTTdetId));

    if ((ttFound != theTTlist_.end()) || *(theTTlist_.end()) == thisTTdetId) {
      output.push_back(eedf.id());
      EEDataFrame eedf2(output.back());
      std::copy(eedf.frame().begin(), eedf.frame().end(), eedf2.frame().begin());
    }
  }
}

void EcalTBReadout::performReadout(edm::Event &event,
                                   const EcalTrigTowerConstituentsMap &theTTmap,
                                   const EBDigiCollection &input,
                                   EBDigiCollection &output) {
  // TB readout
  // step 1: get the target crystal index

  edm::Handle<PEcalTBInfo> theEcalTBInfo;
  event.getByLabel(ecalTBInfoLabel_, theEcalTBInfo);

  int crysId = theEcalTBInfo->nCrystal();

  // step 2: update (if needed) the TT list to be read

  findTTlist(crysId, theTTmap);

  // step 3: perform the readout

  readOut(input, output, theTTmap);
}

void EcalTBReadout::performReadout(edm::Event &event,
                                   const EcalTrigTowerConstituentsMap &theTTmap,
                                   const EEDigiCollection &input,
                                   EEDigiCollection &output) {
  // TB readout
  // step 1: get the target crystal index

  edm::Handle<PEcalTBInfo> theEcalTBInfo;
  event.getByLabel(ecalTBInfoLabel_, theEcalTBInfo);

  int crysId = theEcalTBInfo->nCrystal();

  // step 2: update (if needed) the TT list to be read

  findTTlist(crysId, theTTmap);

  // step 3: perform the readout

  readOut(input, output, theTTmap);
}