LArCalibPatchingAlg.icc

Go to the documentation of this file.

//Dear emacs, this is -*-c++-*-
/*
  Copyright (C) 2002-2022 CERN for the benefit of the ATLAS collaboration
*/
 
#include <tuple>
 
template<class CONDITIONSCONTAINER>
LArCalibPatchingAlg<CONDITIONSCONTAINER>::LArCalibPatchingAlg (const std::string& name, ISvcLocator* pSvcLocator) : 
  AthAlgorithm(name,pSvcLocator),
  m_onlineHelper(0),
  m_caloId(0),
  m_contIn(0),
  m_contOut(0),
  m_patchMethod(PhiAverage) { }
 
 
template<class CONDITIONSCONTAINER>
StatusCode LArCalibPatchingAlg<CONDITIONSCONTAINER>::initialize() {
 
  ATH_CHECK( m_BCKey.initialize() );
  ATH_CHECK( m_cablingKey.initialize() );
  ATH_CHECK( m_CLKey.initialize() );
 
 
  if(m_patchMethodProp=="FEBNeighbor") {
    m_patchMethod=FEBNeighbor;
    return StatusCode::SUCCESS;
  }
  else if (m_patchMethodProp=="PhiNeighbor") {
    m_patchMethod=PhiNeighbor;
    return StatusCode::SUCCESS;
  }
  else if (m_patchMethodProp=="PhiAverage") {
    m_patchMethod=PhiAverage;
    return StatusCode::SUCCESS;
  }
  else if (m_patchMethodProp=="FEBAverage") {
    if (typeid(CONDITIONSCONTAINER)==typeid(LArCaliWaveContainer)) {
      ATH_MSG_ERROR ( "FEBAverage not implemented for CaliWaveContainer." 
                      << "Please choose other patching strategy" );
      return StatusCode::FAILURE;
    }
    m_patchMethod=FEBAverage;
    return StatusCode::SUCCESS;
  } else if (m_patchMethodProp=="SetZero") {
    if (typeid(CONDITIONSCONTAINER)==typeid(LArCaliWaveContainer)) {
      ATH_MSG_ERROR ( "SetZero not implemented for CaliWaveContainer." 
                      << "Please choose other patching strategy" );
      return StatusCode::FAILURE;
    }
    m_patchMethod=SetZero;
    return StatusCode::SUCCESS;
  }
 
  ATH_MSG_ERROR ( "Unknown patching method: " << m_patchMethodProp );
  ATH_MSG_ERROR ( "Allowed values: [Empty, FEBNeighbor, PhiNeighbor, PhiAverage, SetZero]" );
  
  return StatusCode::FAILURE; 
}
 
template<class CONDITIONSCONTAINER>
StatusCode LArCalibPatchingAlg<CONDITIONSCONTAINER>::stop() {
  ATH_MSG_INFO ( "Entering LArCalibPatchingAlg" );
  if(m_isSC) {
     const LArOnline_SuperCellID *onlHelper=nullptr;
     ATH_CHECK( detStore()->retrieve(onlHelper, "LArOnline_SuperCellID") );
     m_onlineHelper=dynamic_cast<const LArOnlineID_Base*>(onlHelper);
     const CaloCell_SuperCell_ID *caloid=nullptr;
     ATH_CHECK( detStore()->retrieve(caloid, "CaloCell_SuperCell_ID") );
     m_caloId = dynamic_cast<const CaloCell_Base_ID*>(caloid); 
  } else {
     const LArOnlineID *onlHelper=nullptr;
     ATH_CHECK( detStore()->retrieve(onlHelper, "LArOnlineID") );
     m_onlineHelper=dynamic_cast<const LArOnlineID_Base*>(onlHelper);
     const CaloCell_ID *caloid=nullptr;
     ATH_CHECK( detStore()->retrieve(caloid, "CaloCell_ID") );
     m_caloId = dynamic_cast<const CaloCell_Base_ID*>(caloid); 
  }
 
  if(m_isSC) m_bcMask.setSC();
  ATH_CHECK(m_bcMask.buildBitMask(m_problemsToPatch,msg()));
 
  const EventContext& ctx = Gaudi::Hive::currentContext();
 
  SG::ReadCondHandle<LArBadChannelCont> readHandle{m_BCKey, ctx};
  const LArBadChannelCont *bcCont {*readHandle};
  if(!bcCont) {
     ATH_MSG_ERROR( "Do not have Bad chan container !!!" );
     return StatusCode::FAILURE;
  }
 
  SG::ReadCondHandle<LArOnOffIdMapping> cablingHdl{m_cablingKey, ctx};
  const LArOnOffIdMapping* cabling = *cablingHdl;
  if(!cabling) {
     ATH_MSG_ERROR( "Do not have OnOff Id mapping !!!" );
     return StatusCode::FAILURE;
  }
 
  SG::ReadCondHandle<LArCalibLineMapping> clHdl{m_CLKey, ctx};
  const LArCalibLineMapping *clCont {*clHdl};
  if(!clCont) {
     ATH_MSG_ERROR( "Do not have calib line mapping !!!" );
     return StatusCode::FAILURE;
  }
 
  if (m_newContainerKey.size()) {
    //New container key give -> different containers for reading and writing
    ATH_CHECK( detStore()->retrieve(m_contIn,m_containerKey) ); //const-retrieve
    
    m_contOut=new CONDITIONSCONTAINER();
    m_contOut->setGroupingType((LArConditionsContainerBase::GroupingType)m_contIn->groupingType());
    ATH_CHECK( m_contOut->initialize() );
    ATH_CHECK( detStore()->record(m_contOut,m_newContainerKey) );
    ATH_CHECK( setSymlink(m_contOut) );
    ATH_MSG_INFO ( "Loaded input container " << m_containerKey 
                   << ", write to new container " << m_newContainerKey );
  }
  else { //Same container for reading and writing
    if (m_unlock) {
      ATH_CHECK( detStore()->retrieve(m_contIn,m_containerKey) ); //const-retrieve
      m_contOut=const_cast<CONDITIONSCONTAINER*>(m_contIn);
    }
    else{
      ATH_CHECK( detStore()->retrieve(m_contOut,m_containerKey) ); //non-const retrieve
      m_contIn=const_cast<const CONDITIONSCONTAINER*>(m_contOut);
      ATH_MSG_INFO ( "Work on container '" <<  m_containerKey  << "'" );
    }
  }
 
  LArBadChanBitPacking packing;
  LArBadChanSCBitPacking scpacking;
 
  unsigned maxgain;
  if(m_isSC) maxgain=CaloGain::LARMEDIUMGAIN; else maxgain=CaloGain::LARNGAIN;
  for (unsigned igain=CaloGain::LARHIGHGAIN; 
       igain<maxgain ; ++igain ) {
    CONTIT it=m_contIn->begin(igain);
    CONTIT it_e=m_contIn->end(igain);
    for (;it!=it_e;it++) {
      const HWIdentifier chid = it.channelId();
      if (!cabling->isOnlineConnected(chid)) continue; //Don't care about disconnected channels
      if (m_bcMask.cellShouldBeMasked(bcCont,chid)) {
    const std::string bcType = m_isSC ? scpacking.stringStatus(bcCont->status(chid)) : packing.stringStatus(bcCont->status(chid));
    ATH_MSG_INFO ( "Found problematic channel 0x" << MSG::hex << chid.get_identifier32().get_compact() << MSG::dec << " [" << bcType << "]" 
                       <<" Gain:" << igain << " " <<  m_onlineHelper->channel_name(chid) << ". Trying to patch." );
    if (patch(chid,igain, bcCont, cabling, clCont)) {
      ATH_MSG_INFO ( "Sucessfully patched  channel 0x" << MSG::hex << chid.get_identifier32().get_compact() << MSG::dec <<" Gain:" << igain );
        } else {
      ATH_MSG_WARNING ( "Failed to patch channel 0x" << MSG::hex << chid.get_identifier32().get_compact() << MSG::dec <<" Gain:" << igain );
        } 
      }//end if channel is in bad-channel database
      else 
    if (it->isEmpty()) { //check if data-object is empty (eg the default instance
      ATH_MSG_ERROR ( "The channel 0x" << MSG::hex << chid.get_identifier32().get_compact() << MSG::dec 
           <<" Gain:" << igain << " " << m_onlineHelper->channel_name(chid) 
                          << " has no calibration but is not (yet?) flagged in the bad-channel database" );
      if (m_patchAllMissing) {
        ATH_MSG_INFO ( "Will try to patch anyway." );
        if (patch(chid,igain, bcCont, cabling, clCont)) {
          ATH_MSG_INFO ( "Sucessfully patched  channel 0x" << MSG::hex << chid.get_identifier32().get_compact() << MSG::dec <<" Gain:" << igain );
        } else {
          ATH_MSG_WARNING ( "Failed to patch channel 0x" << MSG::hex << chid.get_identifier32().get_compact() << MSG::dec <<" Gain:" << igain );
        }
      }//end if m_patchAllMissing
      else 
        ATH_MSG_ERROR ( "Channel remains un-patched!" );
    }//end if isEmpty
    }//end loop over all channels
  }//end loop over all gains
  
  std::vector<unsigned> completedChans = m_contOut->completeCorrectionChannels();
  if (completedChans.size()>0 && m_useCorrChannel) {
    msg() << MSG::INFO << "Artificially inserted correction subsets in COOL channels";
    for(size_t j=0;j<completedChans.size();++j) 
      msg() << MSG::INFO << " " << completedChans[j];
    msg() << MSG::INFO << endmsg;
  }
 
 
  ATH_MSG_INFO ( "Done with LArCalibPatchingAlg" );
  ATH_MSG_DEBUG ( detStore()->dump() );
  return StatusCode::SUCCESS;
}
 
template<class CONDITIONSCONTAINER>
bool LArCalibPatchingAlg<CONDITIONSCONTAINER>::patch(const HWIdentifier chid, const int gain, const LArBadChannelCont* bcCont, const LArOnOffIdMapping* cabling, const LArCalibLineMapping *clCont) {
 
  if(m_doNotPatchCBs.size()>0) { // check if this channel is not on excluded CalibBoard
     const std::vector<HWIdentifier>& cLids=clCont->calibSlotLine(chid);
     for(unsigned cl=0; cl<cLids.size(); ++cl) {
        const HWIdentifier calibModuleID = m_onlineHelper->calib_module_Id(cLids[cl]);
        if (std::find(m_doNotPatchCBs.begin(),m_doNotPatchCBs.end(),calibModuleID.get_identifier32().get_compact()) != m_doNotPatchCBs.end()) { // we should not patch this channel
          return false;
        }
     } // over CLs
  }
   
  if (m_patchMethod==FEBNeighbor){
    const HWIdentifier febId=m_onlineHelper->feb_Id(chid);
    const int febChan=m_onlineHelper->channel(chid);
    HWIdentifier chid_patch(0);
    if (febChan>0) {
      //try lower channel
      chid_patch=m_onlineHelper->channel_Id(febId, febChan-1);
      if (cabling->isOnlineConnected(chid_patch) && bcCont->status(chid_patch).good()) {
    const LArCondObj patch=m_contIn->get(chid_patch,gain); //Should be the const-get method
    if (!patch.isEmpty()) {
      StatusCode sc=m_contOut->insertCorrection(chid,patch,gain,m_useCorrChannel); 
      if (sc.isFailure()) 
        ATH_MSG_ERROR ( "Failed to insert correction for channel 0x" << MSG::hex << chid.get_compact() 
                            << MSG::dec <<", gain " << gain << "." );
      else
        ATH_MSG_INFO ( "Replaced channel 0x" << MSG::hex << chid.get_compact() 
                           <<  " by it's left FEB neighbor 0x" << chid_patch.get_compact() << MSG::dec );
      return true;
    }
      }//end if patch connected and good
    }//end if febChan
    if (febChan<(m_onlineHelper->channelInSlotMax(febId)-1)) {
      chid_patch=m_onlineHelper->channel_Id(febId, febChan+1);
      if (cabling->isOnlineConnected(chid_patch) && bcCont->status(chid_patch).good()) {
    const LArCondObj patch=m_contIn->get(chid_patch,gain); //Should be the const-get method
    if (!patch.isEmpty()) {
      StatusCode sc=m_contOut->insertCorrection(chid,patch,gain,m_useCorrChannel); 
      if (sc.isFailure()) {
        ATH_MSG_ERROR ( "Failed to insert correction for channel 0x" << MSG::hex << chid.get_compact() 
                            << MSG::dec << ", gain " << gain << "." );
      }
      else {
        ATH_MSG_INFO ( "Replaced channel 0x" << MSG::hex << chid.get_compact() 
                           <<  " by it's right FEB neighbor 0x" << chid_patch.get_compact() << MSG::dec );
        return true;
      }
    }//end if patch connected and good
      }//end if chan<max
    }
    ATH_MSG_ERROR ( "None of the FEB neighbors is good!" );
    return false;
  } else if (m_patchMethod==PhiNeighbor) {
    //    (*m_log) << MSG::ERROR << "Patching Method 'Phi-neighbor' not yet implemented." << endmsg;
    try {
      const Identifier id=cabling->cnvToIdentifier(chid);
      int eta, phi, phi_min, phi_max, phi_range;
      Identifier regionID;
      regionID=m_caloId->region_id(id);
      eta      = m_caloId->eta(id);
      phi      = m_caloId->phi(id);
      phi_min =m_caloId->phi_min(regionID);
      phi_max =m_caloId->phi_max(regionID);
      phi_range=phi_max-phi_min+1;
      if (eta==CaloCell_ID::NOT_VALID || phi==CaloCell_ID::NOT_VALID ||  phi_min==CaloCell_ID::NOT_VALID ||  phi_max==CaloCell_ID::NOT_VALID) {
            ATH_MSG_ERROR ( "CaloCell_ID returned NOT_VALID for offline id 0x"<< std::hex << id.get_compact() 
                          <<", online id 0x" << chid.get_compact() << std::dec );
            return false;
      }
      ATH_MSG_VERBOSE ( "Problem channel has phi="<< phi << " eta=" << eta );
      //Try both phi-neighbors
      int phi_list[4]={phi-1,phi+1,phi-2,phi+2};
      for (unsigned i=0;i<4;i++) {
    int phi_patch=phi_list[i];
    if (phi_patch<m_caloId->phi_min(regionID)) phi_patch=phi_patch+phi_range;
    if (phi_patch>m_caloId->phi_max(regionID)) phi_patch=phi_patch-phi_range+phi_min;
    ATH_MSG_VERBOSE ( "Iteration " << i << " Using cell with phi=" 
                          << phi_patch << " eta=" << eta );
    //std::cout << "i=" << i << " Using cell with phi=" 
    //    << phi_patch << " eta=" << eta << std::endl;
    Identifier patch_id=m_caloId->cell_id(regionID,eta,phi_patch);
    HWIdentifier chid_patch=cabling->createSignalChannelID(patch_id);
    if (bcCont->status(chid_patch).good()) {
      const LArCondObj patch=m_contIn->get(chid_patch,gain); 
      if (!patch.isEmpty()) {
        StatusCode sc=m_contOut->insertCorrection(chid,patch,gain,m_useCorrChannel); 
        if (sc.isFailure())
          ATH_MSG_ERROR ( "Failed to insert correction for channel 0x" << MSG::hex << chid.get_compact() 
                              << MSG::dec << ", gain " << gain << "." );
         else
           ATH_MSG_INFO ( "Replaced channel 0x" << MSG::hex << chid.get_compact() 
                              <<  " by neighbor 0x" << chid_patch.get_compact() << MSG::dec << " (phi " << phi << " to " << phi_patch <<")" );
         return true;
       }//end if isEmtpy()
     }//end if neighbor is good
      }// end for loop
      ATH_MSG_ERROR ( "All phi-neighbors of channel 0x" << MSG::hex << chid.get_compact() << MSG::dec 
                      << " are either absent or bad." );
      return false; 
    }catch(LArID_Exception& except) {
      ATH_MSG_ERROR ( "LArID_Exception caught!" );
      return false;
    }
  }  
  else if (m_patchMethod==PhiAverage) {
    LArCondObj patch;
    if (!getAverage(chid,gain,patch,bcCont,cabling)){
      ATH_MSG_ERROR ( "Failed get phi-average!" );
      return false;
    } else {
      ATH_MSG_DEBUG ( "Got a phi-average..." );
    }
    ATH_CHECK( m_contOut->insertCorrection(chid,patch,gain,m_useCorrChannel), false );
    return true;
  }
  else if (m_patchMethod==FEBAverage) {
    LArCondObj patch;
    if (!getAverage(chid,gain,patch,bcCont, cabling, false)){
      ATH_MSG_ERROR ( "Failed get FEB-average!" );
      return false;
    } else {
      ATH_MSG_DEBUG ( "Got a FEB-average..." );
    }
    ATH_CHECK( m_contOut->insertCorrection(chid,patch,gain,m_useCorrChannel), false );
    return true;
  }
  else if (m_patchMethod==SetZero) {
    LArCondObj patch;
    if (!setZero(chid, gain, patch)){
      ATH_MSG_ERROR ( "Failed set Zero!" );
      return false;
    } else {
      ATH_MSG_DEBUG ( "Set zero ..." );
    }
    ATH_CHECK( m_contOut->insertCorrection(chid,patch,gain,m_useCorrChannel), false );
    return true;
  }
  else //failed...
    ATH_MSG_ERROR ( "Unknown correction method." );
  return false;
}
 
 
template<class CONDITIONSCONTAINER>
std::vector<HWIdentifier>& LArCalibPatchingAlg<CONDITIONSCONTAINER>::getPhiRing(const HWIdentifier chid, const LArBadChannelCont* bcCont, const LArOnOffIdMapping* cabling, unsigned distance) {
  if (distance==0) distance=1;
  m_idList.clear();
  try {
    const Identifier id=cabling->cnvToIdentifier(chid);
    int eta,phi, phi_min, phi_max;
    Identifier regionID;
    regionID=m_caloId->region_id(id);
    eta      = m_caloId->eta(id);
    phi      = m_caloId->phi(id);
    phi_min =m_caloId->phi_min(regionID);
    phi_max =m_caloId->phi_max(regionID);
    ATH_MSG_VERBOSE ( "Assembling phi-ring for eta=" << eta << " phi=" << phi );
    //std::cout << "Assembling phi-ring for eta=" << eta << " phi=" << phi << std::endl;
    if (eta==CaloCell_ID::NOT_VALID || phi==CaloCell_ID::NOT_VALID ||  phi_min==CaloCell_ID::NOT_VALID ||  phi_max==CaloCell_ID::NOT_VALID) {
      ATH_MSG_ERROR ( "CaloCell_ID returned NOT_VALID for offline id 0x"<< std::hex << id.get_compact() 
                      <<", online id 0x" << chid.get_compact() << std::dec );
      return m_idList;
    }
    
    if ((phi_max-phi_min)%distance) {
      ATH_MSG_ERROR ( "Can't divide " << (phi_min-phi_max) << " by " << distance );
      return m_idList;
    }
 
    int nSteps=(phi_max-phi_min)/distance;
    
    for (int i=1;i<=nSteps;i++) {
      int phi_patch=phi+i*distance;
      if (phi_patch>phi_max) phi_patch=phi_patch-phi_max+phi_min-1;
      ATH_MSG_VERBOSE ( "i=" << i << " Adding cell with phi=" 
                        << phi_patch << " eta=" << eta );
      //std::cout << "i=" << i << " Adding cell with phi=" 
      //    << phi_patch << " eta=" << eta << std::endl;
      Identifier patch_id=m_caloId->cell_id(regionID,eta,phi_patch);
      HWIdentifier chid_patch=cabling->createSignalChannelID(patch_id);
      if (bcCont->status(chid_patch).good()) {
    m_idList.push_back(chid_patch);
      }
      else
        ATH_MSG_VERBOSE ( "This cell is bad as well. Ignored." );
    }//end loop over phi-steps
  
  }catch(LArID_Exception& except) {
    ATH_MSG_ERROR ( "LArID_Exception caught!" );
  }  
 return m_idList;
}
 
template<class CONDITIONSCONTAINER>
std::vector<HWIdentifier>& LArCalibPatchingAlg<CONDITIONSCONTAINER>::getFEBChans(const HWIdentifier chid, const LArBadChannelCont* bcCont) {
  m_idList.clear();
  HWIdentifier febid = m_onlineHelper->feb_Id(chid);
  ATH_MSG_VERBOSE ( "Assembling list of channels for FEB=" << febid  );
  for (int i=0;i<128;++i) {
      HWIdentifier fchan = m_onlineHelper->channel_Id(febid,i);
      if(fchan == chid) continue;
      ATH_MSG_VERBOSE  ( " Adding channel =" << i );
      if (bcCont->status(fchan).good()) {
    m_idList.push_back(fchan);
      } else {
        ATH_MSG_VERBOSE ( "This channel is bad as well. Ignored." );
      }
  }//end loop over chans
  
  return m_idList;
}
 
 
template<class CONDITIONSCONTAINER>
bool LArCalibPatchingAlg<CONDITIONSCONTAINER>::getAverage(const HWIdentifier chid, const int gain, LArRampP1& patch, const LArBadChannelCont* bcCont, const LArOnOffIdMapping* cabling, bool isphi) {
 
  std::vector<HWIdentifier>& symCells =
    isphi ? getPhiRing(chid, bcCont, cabling) : getFEBChans(chid, bcCont);
  if (symCells.empty()) {
    ATH_MSG_ERROR ( "No symmetry cells found!" );
    return false;
  }
  size_t s=m_contIn->get(symCells[0],gain).m_vRamp.size();
  patch.m_vRamp.clear();
  patch.m_vRamp.resize(s);
  unsigned nCells=0;
  for (HWIdentifier hwid : symCells) {
    const LArRampP1& ramp=m_contIn->get(hwid,gain);
    if (ramp.m_vRamp.size()==0) continue; //This one is empty...
    if (ramp.m_vRamp.size()!=s) {
      if(isphi) {
        ATH_MSG_WARNING ("Cell with same phi but different size of ramp polynom found!" );
      } else {
        ATH_MSG_WARNING ("Cell with same FEB but different size of ramp polynom found!" );
      }
      continue;
    }
    msg() << MSG::DEBUG << "Adding cell 0x"<< std::hex << hwid.get_compact() << std::dec << " Ramp:";
    for (size_t i=0;i<s;i++) {
      patch.m_vRamp[i]+=ramp.m_vRamp[i];
      msg() << MSG::DEBUG << ramp.m_vRamp[i] << " ";
    }
    msg() << MSG::DEBUG << endmsg;
    nCells++;
  }
  if (nCells==0) {
    if(isphi) {
      ATH_MSG_ERROR ( "No good ramp with same phi found!" );
    } else {
      ATH_MSG_ERROR ( "No good ramp with same FEB found!" );
    }
    return false;
  }
 for (size_t i=0;i<s;i++)
   patch.m_vRamp[i]=patch.m_vRamp[i]/nCells;
 
 //FIXME: We should somehow watch the rms....
 msg() << MSG::INFO << "Patched Ramp (based on " << nCells << " channels):" ; 
 for (size_t i=0;i<s;i++)
   msg() << MSG::INFO << " " << patch.m_vRamp[i];
 msg() << MSG::INFO << endmsg;
 return true;
}
 
 
template<class CONDITIONSCONTAINER>
bool LArCalibPatchingAlg<CONDITIONSCONTAINER>::getAverage(const HWIdentifier chid, const int gain, LArOFCP1& patch, const LArBadChannelCont* bcCont, const LArOnOffIdMapping* cabling, bool /*isphi*/) {
 
  std::vector<HWIdentifier>& symCells=getPhiRing(chid, bcCont, cabling);  
 
  if (symCells.empty()) {
    ATH_MSG_ERROR ( "No symmetry cells found!" );
    return false;
  }
 
  const size_t nPhases=m_contIn->get(symCells[0],gain).m_vOFC_a.size();
  if (!nPhases) {
    ATH_MSG_ERROR ( "OFC of neighbor nPhase=0!" );
    return false;
  }
  const size_t nSamples=m_contIn->get(symCells[0],gain).mvOFC_a[0].size();
  if (!nSamples) {
    ATH_MSG_ERROR ( "OFC of neighbor nSamples=0!" );
    return false;
  }
 
  float timeOffset = m_contIn->get(symCells[0],gain).m_timeOffset;
  float timeBinWidth = m_contIn->get(symCells[0],gain).m_timeBinWidth;
  std::vector<std::vector<float> > ofc_a;
  std::vector<std::vector<float> > ofc_b;
  ofc_a.resize(nPhases,std::vector<float>(nSamples));
  ofc_b.resize(nPhases,std::vector<float>(nSamples));
 
  unsigned nCells=0;
  for (HWIdentifier hwid : symCells) {
    const LArOFCP1& ofc=m_contIn->get(hwid,gain);
    if (ofc.OFC_aSize()==0 || ofc.OFC_bSize()==0) continue; //This one is empty...
    if (ofc.OFC_aSize()!=nPhases) {
      ATH_MSG_WARNING ("Cell with same phi but different nPhases found! Ignored" );
      continue;
    }
    if (ofc.timeOffset()!=timeOffset) {
      ATH_MSG_WARNING ("Cell with same phi but different time-offset found! Ignored" );
      continue;
    }
    if (ofc.timeBinWidth()!=timeBinWidth) {
      ATH_MSG_WARNING ("Cell with same phi but different time-offset found! Ignored" );
      continue;
    }
 
    ATH_MSG_DEBUG ( "Adding cell 0x"<< std::hex << hwid.get_compact() << std::dec );
    for (size_t iPhase=0;iPhase<nPhases;++iPhase) {
      //Check size of vector?
      for (size_t iSample=0;iSample<nSamples;++iSample) {
    ofc_a[iPhase][iSample]+=ofc.OFC_a(iPhase)[iSample];
    ofc_b[iPhase][iSample]+=ofc.OFC_b(iPhase)[iSample];
      }
    }
    nCells++;
  }
  if (nCells==0) {
    ATH_MSG_ERROR ( "No good OFC set with same phi found!" );
    return false;
  }
 
  for (size_t iPhase=0;iPhase<nPhases;++iPhase) {
    //Check size of vector?
    for (size_t iSample=0;iSample<nSamples;++iSample) {
      ofc_a[iPhase][iSample]/=nCells;
      ofc_b[iPhase][iSample]/=nCells;
    }
  }
 //FIXME: We should somehow watch the rms....
 
  LArOFCP1 tmp (timeOffset,
                timeBinWidth,
                ofc_a, ofc_b);
  patch.setFrom (tmp);
 
  return true;
}
 
template<class CONDITIONSCONTAINER>
bool LArCalibPatchingAlg<CONDITIONSCONTAINER>::getAverage(const HWIdentifier chid, const int gain, LArCaliWaveVec& patch,const LArBadChannelCont* bcCont, const LArOnOffIdMapping* cabling, bool /*isphi*/) {
 
  const std::vector<HWIdentifier>& symCells=getPhiRing(chid, bcCont, cabling);  
  
 
  if (symCells.empty()) {
    ATH_MSG_ERROR ( "No symmetry cells found!" );
    return false;
  }
  patch.clear();
 
  struct perDAC_t {
    std::vector<std::tuple<HWIdentifier,double,const LArCaliWave*> > tmax_wave; //Identifier, tmax and ptr to wave, all neighbors
    double tmaxAvg=0; //average peak time of all neighbors
  };
 
  LArWaveHelper wHelper;
 
  // Organize waves per DAC and calculate tmax (peaking time)
  std::map<int,perDAC_t> neighbors_per_dac; //In most cases, there is only one DAC, eg one entry in this map
 
  ATH_MSG_DEBUG(" symCells size "<<symCells.size());
  for (HWIdentifier hwid : symCells) {
    for (const LArCaliWave& cwave : m_contIn->get(hwid,gain)) {
      ATH_MSG_VERBOSE("found wave in channel " << m_onlineHelper->channel_name(hwid) << ", gain " << gain);
      const int ourDAC = cwave.getDAC();
      double tmax=-1;
      wHelper.getDMax(cwave, tmax);
      if (tmax > 0 ) {  
    auto& dacWave=neighbors_per_dac[ourDAC];
    dacWave.tmax_wave.emplace_back(hwid,tmax,&cwave);
    dacWave.tmaxAvg+=tmax;
      }
      else {
    ATH_MSG_WARNING("Ignoring wave with peak-time=" << tmax << ", DAC=" << ourDAC << " found in channel " << m_onlineHelper->channel_name(hwid) << ", gain " << gain);
      }
    }//end loop over DACs
  }// end loop over sym-cells
  
  ATH_MSG_DEBUG(" neighbors_per_dac size "<<neighbors_per_dac.size());
 
  for (auto& dacWave : neighbors_per_dac) {
    dacWave.second.tmaxAvg/=dacWave.second.tmax_wave.size();
    ATH_MSG_DEBUG("Average tmax computed " << dacWave.second.tmaxAvg << " for DAC=" << dacWave.first << ", to patch channel " << m_onlineHelper->channel_name(chid) << ", gain " << gain);
  }
  
  //Align waves accoring to their peak-time and sum them
  std::vector<std::pair<LArCaliWave,unsigned> > alignedWaveSums;
  
  for (const auto& dacWave : neighbors_per_dac) {//loop over DAC values
    const int& dac=dacWave.first;
    const perDAC_t& neighbors=dacWave.second;
    const double& tmaxAvg=neighbors.tmaxAvg;
    bool first=true;
 
    for (const auto& [hwid, tmax, wave] : neighbors.tmax_wave) {
      if (std::fabs(tmaxAvg-tmax) > wave->getWave().size()*wave->getDt()/2.0) {//Don't allow very long time shifts  
    ATH_MSG_WARNING("Peaking time of symmetric channel " << m_onlineHelper->channel_name(hwid) 
            << ", DAC=" << dac << ", gain=" << gain << " is too far off the average (" << std::fabs(tmaxAvg-tmax) << "). Ignoring.");
    continue;
      }
      if (first) {
    //First (useable) wave for this DAC value. Initialize a LArCaliWave object for it
    alignedWaveSums.emplace_back(LArCaliWave(wHelper.Dtranslate(*wave, tmaxAvg - tmax).getWave(),wave->getDt(),wave->getDAC(),wave->getIsPulsedInt(),wave->getFlag()),1);
    ATH_MSG_DEBUG("Adding wave of symmetric channel " << m_onlineHelper->channel_name(hwid) 
            << ", DAC=" << dac << ", gain=" << gain << " shifted by " << tmaxAvg-tmax);
    first=false;
      }
      else {
    //Subsequent waves fro this DAC value: Add them. 
    alignedWaveSums.back().first+=wHelper.Dtranslate(*wave, tmaxAvg - tmax);
    alignedWaveSums.back().second++;
    ATH_MSG_DEBUG("Summing wave of symmetric channel " << m_onlineHelper->channel_name(hwid) 
              << ", DAC=" << dac << ", gain=" << gain << " shifted by " << tmaxAvg-tmax);
      }
    }//end loop over symmetic waves for one dac value
    if (first) {
      ATH_MSG_ERROR("Cannot calculate patch wave for channel " <<  m_onlineHelper->channel_name(chid)  << ", gain="<< gain << ", DAC=" << dac << ", no usable neighbors found");
      return false;
    }
    
  }//end loop over DAC values
 
  //Final loop to divide wave by number-of-waves to get the average, then push-back into return container
  for (auto& waveSum : alignedWaveSums) { //Loop over DAC values
    LArCaliWave& wave=waveSum.first;
    wave*=(1./waveSum.second);
    ATH_MSG_INFO("Calculated average wave for channel " <<  m_onlineHelper->channel_name(chid)  << ", gain="<< gain << ", DAC=" << waveSum.first.getDAC() 
         << ", based on " << waveSum.second << " neighbors");
    patch.push_back(wave);
  }
  return true;
}
 
#ifdef LARRAWCONDITIONS_LARMPHYSOVERMCALP
template<class CONDITIONSCONTAINER>
bool LArCalibPatchingAlg<CONDITIONSCONTAINER>::getAverage(const HWIdentifier chid, const int gain, LArMphysOverMcalP1& patch, const LArBadChannelCont* bcCont, const LArOnOffIdMapping* cabling, bool isphi) {
 
  std::vector<HWIdentifier>& symCells=getPhiRing(chid, bcCont, cabling);  
 
  patch.m_MphysOverMcal=0;
 
  unsigned nCells=0;
  for (HWIdentifier hwid : symCells) {
    const float mPmC=m_contIn->get(hwid,gain);
    if (mPmC>0) {
      patch.m_MphysOverMcal+=mPmC;
      nCells++;
    }
  }
 
 if (nCells==0) {
   ATH_MSG_ERROR ( "No good symmetry cells found!" );
   return false;
 }
 patch.m_MphysOverMcal/=nCells;
 return true;
}
#endif
 
 
 
#ifdef LARRAWCONDITIONS_LARSINGLEFLOATP
template<class CONDITIONSCONTAINER>
bool LArCalibPatchingAlg<CONDITIONSCONTAINER>::getAverage(const HWIdentifier chid, const int gain, LArSingleFloatP& patch, const LArBadChannelCont* bcCont, const LArOnOffIdMapping* cabling, bool /*isphi*/) {
 
  std::vector<HWIdentifier>& symCells=getPhiRing(chid, bcCont, cabling);  
 
  patch.m_data=0;
 
  unsigned nCells=0;
  for (HWIdentifier hwid : symCells) {
    const LArSingleFloatP& sf=m_contIn->get(hwid,gain);
    if (!sf.isEmpty()) {
      patch.m_data+=sf.m_data;
      nCells++;
    }
  }
 
 if (nCells==0) {
   ATH_MSG_ERROR ( "No good symmetry cells found!" );
   return false;
 }
 patch.m_data/=nCells;
 return true;
}
#endif
 
template<class CONDITIONSCONTAINER>
bool LArCalibPatchingAlg<CONDITIONSCONTAINER>::setZero(const HWIdentifier chid, const int gain, LArRampP1& patch) {
 
  size_t s=m_contIn->get(chid,gain).m_vRamp.size();
  patch.m_vRamp.clear();
  patch.m_vRamp.resize(s);
  for (size_t i=0; i<s; ++i) patch.m_vRamp[i] = 0.;
  return true;
}
 
template<class CONDITIONSCONTAINER>
bool LArCalibPatchingAlg<CONDITIONSCONTAINER>::setZero(const HWIdentifier chid, const int gain, LArOFCP1& patch) {
 
  const size_t nPhases=m_contIn->get(chid,gain).m_vOFC_a.size();
  if (!nPhases) {
    ATH_MSG_ERROR ( "OFC of nPhase=0 !" );
    return false;
  }
  const size_t nSamples=m_contIn->get(chid,gain).mvOFC_a[0].size();
  if (!nSamples) {
    ATH_MSG_ERROR ( "OFC of nSamples=0 !" );
    return false;
  }
 
  LArOFCP1 newvals (patch.timeOffset(), patch.timeBinWidth(),
                    std::vector<std::vector<float> > (nPhases, std::vector<float> (nSamples, 0)),
                    std::vector<std::vector<float> > (nPhases, std::vector<float> (nSamples, 0)));
  patch.setFrom (newvals);
  return true;
}
 
template<class CONDITIONSCONTAINER>
bool LArCalibPatchingAlg<CONDITIONSCONTAINER>::setZero(const HWIdentifier /*chid*/, const int /*gain*/, LArCaliWaveVec& /*patch*/) {
  ATH_MSG_ERROR ( "Not implemented, should not come here !!!");
  return false;
}
 
#ifdef LARRAWCONDITIONS_LARSINGLEFLOATP
template<class CONDITIONSCONTAINER>
bool LArCalibPatchingAlg<CONDITIONSCONTAINER>::setZero(LArSingleFloatP& patch) {
  patch.m_data=0.;
  return true;
}
#endif
 
 
template<class CONDITIONSCONTAINER>
StatusCode  LArCalibPatchingAlg<CONDITIONSCONTAINER>::setSymlink(const LArRampComplete* ramp) const {
  ATH_CHECK( detStore()->symLink(ramp, (ILArRamp*)ramp) );
  return StatusCode::SUCCESS;
}
 
 
template<class CONDITIONSCONTAINER>
StatusCode  LArCalibPatchingAlg<CONDITIONSCONTAINER>::setSymlink(const LArOFCComplete* ofc) const {
  ATH_CHECK( detStore()->symLink(ofc, (ILArOFC*)ofc) );
  return StatusCode::SUCCESS;
}
 
 
template<class CONDITIONSCONTAINER>
StatusCode  LArCalibPatchingAlg<CONDITIONSCONTAINER>::setSymlink(const LArMphysOverMcalComplete* ramp) const {
  ATH_CHECK( detStore()->symLink(ramp, (ILArMphysOverMcal*)ramp) );
  return StatusCode::SUCCESS;
}