Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-09-07 04:37:39

0001 // This is  CSCHitFromStripOnly.cc
0002 
0003 #include "RecoLocalMuon/CSCRecHitD/src/CSCHitFromStripOnly.h"
0004 #include "RecoLocalMuon/CSCRecHitD/src/CSCStripData.h"
0005 #include "RecoLocalMuon/CSCRecHitD/src/CSCStripHitData.h"
0006 #include "RecoLocalMuon/CSCRecHitD/src/CSCStripHit.h"
0007 #include "RecoLocalMuon/CSCRecHitD/src/CSCPedestalChoice.h"
0008 
0009 #include "DataFormats/MuonDetId/interface/CSCDetId.h"
0010 
0011 #include "Geometry/CSCGeometry/interface/CSCLayer.h"
0012 #include "Geometry/CSCGeometry/interface/CSCChamberSpecs.h"
0013 
0014 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0015 #include "FWCore/Utilities/interface/Exception.h"
0016 
0017 #include <algorithm>
0018 #include <string>
0019 #include <vector>
0020 //#include <iostream>
0021 
0022 CSCHitFromStripOnly::CSCHitFromStripOnly(const edm::ParameterSet& ps)
0023     : recoConditions_(nullptr), calcped_(nullptr), ganged_(false) {
0024   useCalib = ps.getParameter<bool>("CSCUseCalibrations");
0025   bool useStaticPedestals = ps.getParameter<bool>("CSCUseStaticPedestals");
0026   int noOfTimeBinsForDynamicPed = ps.getParameter<int>("CSCNoOfTimeBinsForDynamicPedestal");
0027 
0028   theThresholdForAPeak = ps.getParameter<double>("CSCStripPeakThreshold");
0029   theThresholdForCluster = ps.getParameter<double>("CSCStripClusterChargeCut");
0030 
0031   LogTrace("CSCRecHit") << "[CSCHitFromStripOnly] CSCUseStaticPedestals = " << useStaticPedestals;
0032   if (!useStaticPedestals)
0033     LogTrace("CSCRecHit") << "[CSCHitFromStripOnly] CSCNoOfTimeBinsForDynamicPedestal = " << noOfTimeBinsForDynamicPed;
0034 
0035   if (useStaticPedestals) {
0036     calcped_ = new CSCStaticPedestal();
0037   } else {
0038     if (noOfTimeBinsForDynamicPed == 1) {
0039       calcped_ = new CSCDynamicPedestal1();
0040     } else {
0041       calcped_ = new CSCDynamicPedestal2();  // NORMAL DEFAULT!
0042     }
0043   }
0044 }
0045 
0046 CSCHitFromStripOnly::~CSCHitFromStripOnly() { delete calcped_; }
0047 
0048 /* runStrip
0049  *
0050  * Search for strip with ADC output exceeding theThresholdForAPeak.  For each of these strips,
0051  * build a cluster of strip of size theClusterSize (typically 3 strips).  Finally, make
0052  * a Strip Hit out of these clusters by finding the center-of-mass position of the hit
0053  */
0054 std::vector<CSCStripHit> CSCHitFromStripOnly::runStrip(const CSCDetId& id,
0055                                                        const CSCLayer* layer,
0056                                                        const CSCStripDigiCollection::Range& rstripd) {
0057   std::vector<CSCStripHit> hitsInLayer;
0058 
0059   // cache layer info for ease of access
0060   id_ = id;
0061   layer_ = layer;
0062   nstrips_ = layer->chamber()->specs()->nStrips();
0063 
0064   setGanged(false);
0065   if (id_.ring() == 4 && layer_->chamber()->specs()->gangedStrips())
0066     setGanged(true);  //@@ ONLY ME1/1A CAN BE GANGED
0067 
0068   LogTrace("CSCHitFromStripOnly") << "[CSCHitFromStripOnly::runStrip] id= " << id_ << " nstrips= " << nstrips_
0069                                   << " ganged strips? " << ganged();
0070 
0071   tmax_cluster = 5;
0072 
0073   // Get gain correction weights for all strips in layer, and cache in gainWeight.
0074   // They're used in fillPulseHeights below.
0075   // When ME11a is ganged we only need the first 16 values of the 48 filled,
0076   // but 17-48 are just duplicates of 1-16 anyway
0077 
0078   if (useCalib) {
0079     recoConditions_->stripWeights(id, nstrips_, gainWeight);
0080 
0081     // *** START DUMP gainWeight
0082     //    std::cout << "gainWeight for id= " << id_ << " nstrips= " << nstrips_ << std::endl;
0083     //    for ( size_t i = 0; i!=10; ++i ) {
0084     //      for ( size_t j = 0; j!=8; ++j ) {
0085     //        std::cout << gainWeight[i*8 + j] << "   ";
0086     //      }
0087     //      std::cout << std::endl;
0088     //    }
0089     // *** END DUMP gainWeight
0090   }
0091 
0092   // Store pulseheights from SCA and find maxima (potential hits)
0093   fillPulseHeights(rstripd);
0094   findMaxima(id);
0095 
0096   // Make a Strip Hit out of each strip local maximum
0097   for (size_t imax = 0; imax != theMaxima.size(); ++imax) {
0098     // Initialize parameters entering the CSCStripHit
0099     clusterSize = theClusterSize;
0100     theStrips.clear();
0101     strips_adc.clear();
0102     strips_adcRaw.clear();
0103 
0104     // makeCluster calls findHitOnStripPosition to determine the centroid position
0105 
0106     // Remember, the array starts at 0, but the stripId starts at 1...
0107     float strippos = makeCluster(theMaxima[imax] + 1);
0108 
0109     //if ( strippos < 0 || tmax_cluster < 3 ){
0110     // the strippos (as calculated here) is not used later on in
0111     /// fact (20.10.09);
0112     // with the negative charges allowed it can become negative
0113     if (tmax_cluster < 3) {
0114       theClosestMaximum.push_back(99);  // to keep proper vector size
0115       continue;
0116     }
0117     //---- If two maxima are too close the error assigned will be width/sqrt(12) - see CSCXonStrip_MatchGatti.cc
0118     int maximum_to_left = 99;  //---- If there is one maximum - the distance is set to 99 (strips)
0119     int maximum_to_right = 99;
0120     if (imax < theMaxima.size() - 1) {
0121       maximum_to_right = theMaxima.at(imax + 1) - theMaxima.at(imax);
0122     }
0123     if (imax > 0) {
0124       maximum_to_left = theMaxima.at(imax - 1) - theMaxima.at(imax);
0125     }
0126     if (std::abs(maximum_to_right) < std::abs(maximum_to_left)) {
0127       theClosestMaximum.push_back(maximum_to_right);
0128     } else {
0129       theClosestMaximum.push_back(maximum_to_left);
0130     }
0131 
0132     //---- Check if a neighbouring strip is a dead strip
0133     //bool deadStrip = isNearDeadStrip(id, theMaxima.at(imax));
0134     bool deadStripL = isDeadStrip(id, theMaxima.at(imax) - 1, nstrips_);
0135     bool deadStripR = isDeadStrip(id, theMaxima.at(imax) + 1, nstrips_);
0136     short int aDeadStrip = 0;
0137     if (!deadStripL && !deadStripR) {
0138       aDeadStrip = 0;
0139     } else if (deadStripL && deadStripR) {
0140       aDeadStrip = 255;
0141     } else {
0142       if (deadStripL) {
0143         aDeadStrip = theMaxima.at(imax) - 1;
0144       } else {
0145         aDeadStrip = theMaxima.at(imax) + 1;
0146       }
0147     }
0148 
0149     /// L1A (Begin looping)
0150     /// Attempt to redefine theStrips, to encode L1A phase bits
0151     std::vector<int> theL1AStrips;
0152     for (int ila = 0; ila < (int)theStrips.size(); ila++) {
0153       bool stripMatchCounter = false;
0154       for (auto itl1 = rstripd.first; itl1 != rstripd.second; ++itl1) {
0155         int stripNproto = (*itl1).getStrip();
0156         if (!ganged()) {
0157           if (theStrips[ila] == stripNproto) {
0158             stripMatchCounter = true;
0159             auto sz = (*itl1).getOverlappedSample().size();
0160             int L1AbitOnPlace = 0;
0161             for (auto iBit = 0UL; iBit < sz; iBit++) {
0162               L1AbitOnPlace |= ((*itl1).getL1APhase(iBit) << (15 - iBit));
0163             }
0164             theL1AStrips.push_back(theStrips[ila] | L1AbitOnPlace);
0165           }
0166         } else {
0167           for (int tripl = 0; tripl < 3; ++tripl) {
0168             if (theStrips[ila] == (stripNproto + tripl * 16)) {
0169               stripMatchCounter = true;
0170               auto sz = (*itl1).getOverlappedSample().size();
0171               int L1AbitOnPlace = 0;
0172               for (auto iBit = 0UL; iBit < sz; iBit++) {
0173                 L1AbitOnPlace |= ((*itl1).getL1APhase(iBit) << (15 - iBit));
0174               }
0175               theL1AStrips.push_back(theStrips[ila] | L1AbitOnPlace);
0176             }
0177           }
0178         }
0179       }
0180       if (!stripMatchCounter) {
0181         theL1AStrips.push_back(theStrips[ila]);
0182       }
0183     }
0184     /// L1A (end Looping)
0185 
0186     CSCStripHit striphit(id,
0187                          strippos,
0188                          tmax_cluster,
0189                          theL1AStrips,
0190                          strips_adc,
0191                          strips_adcRaw,  /// L1A
0192                          theConsecutiveStrips.at(imax),
0193                          theClosestMaximum.at(imax),
0194                          aDeadStrip);
0195     hitsInLayer.push_back(striphit);
0196   }
0197 
0198   /// Print statement to check StripHit content w/ LA1
0199   /*   
0200       for(std::vector<CSCStripHit>::const_iterator itSHit=hitsInLayer.begin(); itSHit!=hitsInLayer.end(); ++itSHit){
0201          (*itSHit).print(); 
0202          }  
0203   */
0204 
0205   return hitsInLayer;
0206 }
0207 
0208 /* makeCluster
0209  *
0210  */
0211 float CSCHitFromStripOnly::makeCluster(int centerStrip) {
0212   float strippos = -1.;
0213   clusterSize = theClusterSize;
0214   std::vector<CSCStripHitData> stripDataV;
0215 
0216   // We only want to use strip position in terms of strip # for the strip hit. //@@ What other choice is there?
0217 
0218   // If the cluster size is such that you go beyond the edge of detector, shrink cluster appropriately
0219   for (int i = 1; i < theClusterSize / 2 + 1; ++i) {
0220     if (centerStrip - i < 1 || centerStrip + i > int(nstrips_)) {
0221       // Shrink cluster size, but keep it an odd number of strips.
0222       clusterSize = 2 * i - 1;
0223     }
0224   }
0225   for (int i = -clusterSize / 2; i <= clusterSize / 2; ++i) {
0226     CSCStripHitData data = makeStripData(centerStrip, i);
0227     stripDataV.push_back(data);
0228     theStrips.push_back(centerStrip + i);
0229   }
0230   strippos = findHitOnStripPosition(stripDataV, centerStrip);
0231 
0232   LogTrace("CSCHitFromStripOnly") << "[CSCHitFromStripOnly::makeCluster] centerStrip= " << centerStrip
0233                                   << " strippos=" << strippos;
0234 
0235   return strippos;
0236 }
0237 
0238 /** makeStripData
0239  *
0240  */
0241 CSCStripHitData CSCHitFromStripOnly::makeStripData(int centerStrip, int offset) {
0242   CSCStripHitData prelimData;
0243   int thisStrip = centerStrip + offset;
0244 
0245   int tmax = thePulseHeightMap[centerStrip - 1].tmax();
0246   tmax_cluster = tmax;
0247 
0248   std::vector<float> adc(4);
0249   std::vector<float> adcRaw(4);
0250 
0251   // Fill adc & adcRaw
0252 
0253   int istart = tmax - 1;
0254   int istop = std::min(tmax + 2, 7);  // there are only time bins 0-7
0255   adc[3] = 0.1;                       // in case it isn't filled
0256 
0257   if (tmax > 2 && tmax < 7) {  // for time bins 3-6
0258     int ibin = thisStrip - 1;
0259     if (thePulseHeightMap[ibin].valid()) {
0260       std::copy(
0261           thePulseHeightMap[ibin].ph().begin() + istart, thePulseHeightMap[ibin].ph().begin() + istop + 1, adc.begin());
0262 
0263       std::copy(thePulseHeightMap[ibin].phRaw().begin() + istart,
0264                 thePulseHeightMap[ibin].phRaw().begin() + istop + 1,
0265                 adcRaw.begin());
0266     }
0267   } else {
0268     adc[0] = 0.1;
0269     adc[1] = 0.1;
0270     adc[2] = 0.1;
0271     adc[3] = 0.1;
0272     adcRaw = adc;
0273     LogTrace("CSCRecHit") << "[CSCHitFromStripOnly::makeStripData] Tmax out of range: contact CSC expert!";
0274   }
0275 
0276   if (offset == 0) {
0277     prelimData = CSCStripHitData(thisStrip, tmax_cluster, adcRaw, adc);
0278   } else {
0279     int sign = offset > 0 ? 1 : -1;
0280     // If there's another maximum that would like to use part of this cluster,
0281     // it gets shared in proportion to the height of the maxima
0282     for (int i = 1; i <= clusterSize / 2; ++i) {
0283       // Find the direction of the offset
0284       int testStrip = thisStrip + sign * i;
0285       std::vector<int>::iterator otherMax = find(theMaxima.begin(), theMaxima.end(), testStrip - 1);
0286 
0287       // No other maxima found, so just store
0288       if (otherMax == theMaxima.end()) {
0289         prelimData = CSCStripHitData(thisStrip, tmax_cluster, adcRaw, adc);
0290       } else {
0291         // Another maximum found - share
0292         std::vector<float> adc1(4);
0293         std::vector<float> adcRaw1(4);
0294         std::vector<float> adc2(4);
0295         std::vector<float> adcRaw2(4);
0296         // In case we only copy (below) into 3 of the 4 bins i.e. when istart=5, istop=7
0297         adc1[3] = 0.1;
0298         adc2[3] = 0.1;
0299         adcRaw1[3] = 0.1;
0300         adcRaw2[3] = 0.1;
0301 
0302         // Fill adcN with content of time bins tmax-1 to tmax+2 (if it exists!)
0303         if (tmax > 2 && tmax < 7) {  // for time bin tmax from 3-6
0304           int ibin = testStrip - 1;
0305           int jbin = centerStrip - 1;
0306           if (thePulseHeightMap[ibin].valid()) {
0307             std::copy(thePulseHeightMap[ibin].ph().begin() + istart,
0308                       thePulseHeightMap[ibin].ph().begin() + istop + 1,
0309                       adc1.begin());
0310             std::copy(thePulseHeightMap[ibin].phRaw().begin() + istart,
0311                       thePulseHeightMap[ibin].phRaw().begin() + istop + 1,
0312                       adcRaw1.begin());
0313           }
0314 
0315           if (thePulseHeightMap[jbin].valid()) {
0316             std::copy(thePulseHeightMap[jbin].ph().begin() + istart,
0317                       thePulseHeightMap[jbin].ph().begin() + istop + 1,
0318                       adc2.begin());
0319 
0320             std::copy(thePulseHeightMap[jbin].phRaw().begin() + istart,
0321                       thePulseHeightMap[jbin].phRaw().begin() + istop + 1,
0322                       adcRaw2.begin());
0323           }
0324         } else {
0325           adc1.assign(4, 0.1);
0326           adcRaw1 = adc1;
0327           adc2.assign(4, 0.1);
0328           adcRaw2 = adc2;
0329         }
0330 
0331         // Scale shared strip B ('adc') by ratio of peak of ADC counts from central strip A ('adc2')
0332         // to sum of A and neighbouring maxima C ('adc1')
0333 
0334         for (size_t k = 0; k < 4; ++k) {
0335           if (adc1[k] > 0 && adc2[k] > 0)
0336             adc[k] = adc[k] * adc2[k] / (adc1[k] + adc2[k]);
0337           if (adcRaw1[k] > 0 && adcRaw2[k] > 0)
0338             adcRaw[k] = adcRaw[k] * adcRaw2[k] / (adcRaw1[k] + adcRaw2[k]);
0339         }
0340         prelimData = CSCStripHitData(thisStrip, tmax_cluster, adcRaw, adc);
0341       }
0342     }
0343   }
0344   return prelimData;
0345 }
0346 
0347 /* fillPulseHeights
0348  *
0349  */
0350 void CSCHitFromStripOnly::fillPulseHeights(const CSCStripDigiCollection::Range& rstripd) {
0351   // Loop over strip digis in one CSCLayer and fill PulseHeightMap with pedestal-subtracted
0352   // SCA pulse heights.
0353 
0354   for (auto& ph : thePulseHeightMap)
0355     ph.reset();
0356 
0357   // for storing sca pulseheights once they may no longer be integer (e.g. after ped subtraction)
0358   for (CSCStripDigiCollection::const_iterator it = rstripd.first; it != rstripd.second; ++it) {
0359     int thisChannel = (*it).getStrip();
0360     auto& stripData = thePulseHeightMap[thisChannel - 1];
0361     auto& scaRaw = stripData.phRaw_;
0362     auto& sca = stripData.ph_;
0363 
0364     auto const& scaOri = (*it).getADCCounts();
0365     assert(scaOri.size() == 8);
0366     // Fill sca from scaRaw, implicitly converting to float
0367     std::copy(scaOri.begin(), scaOri.end(), scaRaw.begin());
0368     std::copy(scaRaw.begin(), scaRaw.end(), sca.begin());
0369 
0370     //@@ Find bin with largest pulseheight (_before_ ped subtraction - shouldn't matter, right?)
0371     int tmax = std::max_element(sca.begin(), sca.end()) - sca.begin();  // counts from 0
0372 
0373     // get pedestal - calculated as appropriate - for this sca pulse
0374     float ped = calcped_->pedestal(sca, recoConditions_, id_, thisChannel);
0375 
0376     // subtract the pedestal (from BOTH sets of sca pulseheights)
0377     std::for_each(sca.begin(), sca.end(), CSCSubtractPedestal(ped));
0378     std::for_each(scaRaw.begin(), scaRaw.end(), CSCSubtractPedestal(ped));
0379 
0380     //@@ Max in first 3 or last time bins is unacceptable, if so set to zero (why?)
0381     float phmax = 0.f;
0382     if (tmax > 2 && tmax < 7) {
0383       phmax = sca[tmax];
0384     }
0385     stripData.phmax_ = phmax;
0386     stripData.tmax_ = tmax;
0387 
0388     // Fill the map, possibly apply gains from cond data, and unfold ME1A channels
0389     // (To apply gains use CSCStripData::op*= which scales only the non-raw sca ph's;
0390     // but note that both sca & scaRaw are pedestal-subtracted.)
0391 
0392     // From StripDigi, thisChannel labels strip channel. Values phmax, tmax, scaRaw, sca belong to thisChannel
0393     if (useCalib)
0394       stripData *= gainWeight[thisChannel - 1];
0395 
0396     // for ganged ME1a need to duplicate values on istrip=thisChannel to iStrip+16 and iStrip+32
0397     if (ganged()) {
0398       for (int j = 1; j < 3; ++j) {
0399         thePulseHeightMap[thisChannel - 1 + 16 * j] = stripData;
0400       }
0401     }
0402   }
0403 }
0404 
0405 /* findMaxima
0406  *
0407  * fills vector 'theMaxima' with the local maxima in the pulseheight distribution
0408  * of the strips. The threshold defining a maximum is a configurable parameter.
0409  * A typical value is 30.
0410  */
0411 void CSCHitFromStripOnly::findMaxima(const CSCDetId& id) {
0412   theMaxima.clear();
0413   theConsecutiveStrips.clear();
0414   theClosestMaximum.clear();
0415   for (size_t i = 0; i != thePulseHeightMap.size(); ++i) {
0416     // sum 3 strips so that hits between strips are not suppressed
0417     float heightCluster = 0.;
0418 
0419     bool maximumFound = false;
0420     // Left edge of chamber
0421     if (!isDeadStrip(id, i + 1, nstrips_)) {  // Is it i or i+1
0422       if (i == 0) {
0423         heightCluster = thePulseHeightMap[i].phmax() + thePulseHeightMap[i + 1].phmax();
0424         // Have found a strip Hit if...
0425         if (thePulseHeightMap[i].phmax() >= thePulseHeightMap[i + 1].phmax() && isPeakOK(i, heightCluster)) {
0426           maximumFound = true;
0427         }
0428         // Right edge of chamber
0429       } else if (i == thePulseHeightMap.size() - 1) {
0430         heightCluster = thePulseHeightMap[i - 1].phmax() + thePulseHeightMap[i].phmax();
0431         // Have found a strip Hit if...
0432         if (thePulseHeightMap[i].phmax() > thePulseHeightMap[i - 1].phmax() && isPeakOK(i, heightCluster)) {
0433           maximumFound = true;
0434         }
0435         // Any other strips
0436       } else {
0437         heightCluster =
0438             thePulseHeightMap[i - 1].phmax() + thePulseHeightMap[i].phmax() + thePulseHeightMap[i + 1].phmax();
0439         // Have found a strip Hit if...
0440         if (thePulseHeightMap[i].phmax() > thePulseHeightMap[i - 1].phmax() &&
0441             thePulseHeightMap[i].phmax() >= thePulseHeightMap[i + 1].phmax() && isPeakOK(i, heightCluster)) {
0442           maximumFound = true;
0443         }
0444       }
0445     }
0446     //---- Consecutive strips with charge (real cluster); if too wide - measurement is not accurate
0447     if (maximumFound) {
0448       int numberOfConsecutiveStrips = 1;
0449       float testThreshold = 10.;  //---- ADC counts;
0450                                   //---- this is not XTalk corrected so it is correct in first approximation only
0451       int j = 0;
0452       for (int l = 0; l < 8; ++l) {
0453         if (j < 0)
0454           edm::LogWarning("FailedStripCountingWrongConsecutiveStripNumber")
0455               << "This should never occur!!! Contact CSC expert!";
0456         ++j;
0457         bool signalPresent = false;
0458         for (int k = 0; k < 2; ++k) {
0459           j *= -1;  //---- check from left and right
0460           int anotherConsecutiveStrip = i + j;
0461           if (anotherConsecutiveStrip >= 0 && anotherConsecutiveStrip < int(thePulseHeightMap.size())) {
0462             if (thePulseHeightMap[anotherConsecutiveStrip].phmax() > testThreshold) {
0463               ++numberOfConsecutiveStrips;
0464               signalPresent = true;
0465             }
0466           }
0467         }
0468         if (!signalPresent) {
0469           break;
0470         }
0471       }
0472 
0473       bool additional_maxima_found = false;
0474       // search for additional maxima if:
0475       // - hit is closer than 3 strips from the edge
0476       // - enough consecutive strips with signal
0477       // - strip charge distribution looks abnormal
0478 
0479       if (i > 2 && i + 3 < thePulseHeightMap.size() && numberOfConsecutiveStrips > 3) {
0480         //try to look for additional maxima at the left side from the main maxima
0481 
0482         if (((thePulseHeightMap[i + 1].phmax() >= thePulseHeightMap[i - 1].phmax() &&
0483               thePulseHeightMap[i + 1].phmax() >= thePulseHeightMap[i - 2].phmax() &&
0484               thePulseHeightMap[i + 2].phmax() <= thePulseHeightMap[i - 2].phmax()) ||
0485              (thePulseHeightMap[i + 1].phmax() <= thePulseHeightMap[i - 1].phmax() &&
0486               thePulseHeightMap[i + 1].phmax() <= thePulseHeightMap[i - 2].phmax())) &&
0487             //to avoid close maxima delimitation (this is already present in the code)
0488             thePulseHeightMap[i - 1].phmax() >= thePulseHeightMap[i - 2].phmax() &&
0489             //no need in a small charge maxima (might need adjustment)
0490             thePulseHeightMap[i - 2].phmax() > 20) {
0491           additional_maxima_found = true;
0492           theMaxima.push_back(i - 2);  //insert left maxima first
0493           //insert the same number of cosecutive strips, because they belong to both maximas
0494           theConsecutiveStrips.push_back(numberOfConsecutiveStrips);
0495           theMaxima.push_back(i);  //insert main maxima
0496           //insert the same number of cosecutive strips, because they belong to both maximas
0497           theConsecutiveStrips.push_back(numberOfConsecutiveStrips);
0498 
0499         }  //looking for additional maxima on the left
0500 
0501         //try to look for additional maxima at the right side from the main maxima
0502         if (((thePulseHeightMap[i + 1].phmax() >= thePulseHeightMap[i - 1].phmax() &&
0503               thePulseHeightMap[i + 2].phmax() >= thePulseHeightMap[i - 1].phmax()) ||
0504              (thePulseHeightMap[i + 1].phmax() <= thePulseHeightMap[i - 1].phmax() &&
0505               thePulseHeightMap[i + 2].phmax() <= thePulseHeightMap[i - 1].phmax() &&
0506               thePulseHeightMap[i + 2].phmax() >= thePulseHeightMap[i - 2].phmax())) &&
0507             //to avoid close maxima delimitation (this is already present in the code)
0508             thePulseHeightMap[i + 1].phmax() >= thePulseHeightMap[i + 2].phmax() &&
0509             //no need in a small charge maxima (might need adjustment)
0510             thePulseHeightMap[i + 2].phmax() > 20) {
0511           additional_maxima_found = true;
0512           theMaxima.push_back(i);  //insert main maxima first
0513           //insert the same number of cosecutive strips, because they belong to both maximas
0514           theConsecutiveStrips.push_back(numberOfConsecutiveStrips);
0515           theMaxima.push_back(i + 2);  //insert right maxima
0516           //insert the same number of cosecutive strips, because they belong to both maximas
0517           theConsecutiveStrips.push_back(numberOfConsecutiveStrips);
0518 
0519         }  //looking for additional maxima on the right
0520 
0521         //if nothing additional found fill the maxima
0522         if (!additional_maxima_found) {
0523           theMaxima.push_back(i);
0524           theConsecutiveStrips.push_back(numberOfConsecutiveStrips);
0525         }
0526       } else {  //not the case for looking for the additional maxima
0527         theMaxima.push_back(i);
0528         theConsecutiveStrips.push_back(numberOfConsecutiveStrips);
0529       }
0530     }  //if maximafound
0531   }  //all pulses
0532 }  //find maxima procedure
0533 
0534 bool CSCHitFromStripOnly::isPeakOK(int iStrip, float heightCluster) {
0535   int i = iStrip;
0536   bool peakOK = (thePulseHeightMap[i].phmax() > theThresholdForAPeak && heightCluster > theThresholdForCluster &&
0537                  // ... and proper peak time; note that the values below are used elsewhere in this file;
0538                  // they should become parameters or at least constants defined in appropriate place
0539                  thePulseHeightMap[i].tmax() > 2 && thePulseHeightMap[i].tmax() < 7);
0540   return peakOK;
0541 }
0542 
0543 /* findHitOnStripPosition
0544  *
0545  */
0546 float CSCHitFromStripOnly::findHitOnStripPosition(const std::vector<CSCStripHitData>& data, const int& centerStrip) {
0547   float strippos = -1.;
0548 
0549   if (data.empty())
0550     return strippos;
0551 
0552   // biggestStrip is strip with largest pulse height
0553   // Use pointer subtraction
0554 
0555   int biggestStrip = max_element(data.begin(), data.end()) - data.begin();
0556   strippos = data[biggestStrip].strip() * 1.;
0557 
0558   // If more than one strip:  use centroid to find center of cluster
0559   // but only use time bin == tmax (otherwise, bias centroid).
0560   float sum = 0.;
0561   float sum_w = 0.;
0562 
0563   //  std::vector<float> w(4);
0564   // std::vector<float> wRaw(4);
0565 
0566   for (size_t i = 0; i != data.size(); ++i) {
0567     auto const& w = data[i].ph();
0568     auto const& wRaw = data[i].phRaw();
0569 
0570     // (Require ADC to be > 0.)
0571     // No later studies suggest that this only do harm
0572     /*
0573     for ( size_t j = 0; j != w.size(); ++j ) {
0574       if ( w[j] < 0. ) w[j] = 0.001;
0575     }
0576     */
0577 
0578     // Fill the data members
0579     std::copy(w.begin(), w.end(), std::back_inserter(strips_adc));
0580     std::copy(wRaw.begin(), wRaw.end(), std::back_inserter(strips_adcRaw));
0581 
0582     if (data[i].strip() < 1) {
0583       LogTrace("CSCRecHit") << "[CSCHitFromStripOnly::findHitOnStripPosition] problem in indexing of strip, strip= "
0584                             << data[i].strip();
0585     }
0586     sum_w += w[1];
0587     sum += w[1] * data[i].strip();
0588   }
0589 
0590   if (sum_w > 0.)
0591     strippos = sum / sum_w;
0592 
0593   return strippos;
0594 }
0595 
0596 bool CSCHitFromStripOnly::isNearDeadStrip(const CSCDetId& id, int centralStrip, int nstrips) {
0597   //@@ Tim says: not sure I understand this properly... but just moved code to CSCRecoConditions
0598   // where it can handle the conversion from strip to channel etc.
0599   return recoConditions_->nearBadStrip(id, centralStrip, nstrips);
0600 }
0601 
0602 bool CSCHitFromStripOnly::isDeadStrip(const CSCDetId& id, int centralStrip, int nstrips) {
0603   return recoConditions_->badStrip(id, centralStrip, nstrips);
0604 }
0605 
0606 // Define space for static
0607 const int CSCHitFromStripOnly::theClusterSize;