SiTrivialInduceChargeOnStrips Class Reference

#include <SiTrivialInduceChargeOnStrips.h>

Inheritance diagram for SiTrivialInduceChargeOnStrips:

Public Member Functions
void	induce (const SiChargeCollectionDrifter::collection_type &collection_points, const StripGeomDetUnit &det, std::vector< float > &localAmplitudes, size_t &recordMinAffectedStrip, size_t &recordMaxAffectedStrip, const TrackerTopology *tTopo) const
	SiTrivialInduceChargeOnStrips (const edm::ParameterSet &conf, double g)
virtual	~SiTrivialInduceChargeOnStrips ()
Public Member Functions inherited from SiInduceChargeOnStrips
virtual	~SiInduceChargeOnStrips ()

Private Member Functions
void	induceOriginal (const SiChargeCollectionDrifter::collection_type &collection_points, const StripGeomDetUnit &det, std::vector< float > &localAmplitudes, size_t &recordMinAffectedStrip, size_t &recordMaxAffectedStrip, const TrackerTopology *tTopo) const
void	induceVector (const SiChargeCollectionDrifter::collection_type &collection_points, const StripGeomDetUnit &det, std::vector< float > &localAmplitudes, size_t &recordMinAffectedStrip, size_t &recordMaxAffectedStrip, const TrackerTopology *tTopo) const

Private Attributes
const float	geVperElectron
const float	Nsigma
const std::vector< std::vector < float > >	signalCoupling

Detailed Description

Definition at line 9 of file SiTrivialInduceChargeOnStrips.h.

Constructor & Destructor Documentation

SiTrivialInduceChargeOnStrips::SiTrivialInduceChargeOnStrips	(	const edm::ParameterSet &	conf,
		double	g
	)

Definition at line 84 of file SiTrivialInduceChargeOnStrips.cc.

  : signalCoupling(fillSignalCoupling(conf, Ntypes, type)), Nsigma(3.), geVperElectron(g)  {
}

virtual SiTrivialInduceChargeOnStrips::~SiTrivialInduceChargeOnStrips ( )

inlinevirtual

Definition at line 12 of file SiTrivialInduceChargeOnStrips.h.

{}

Member Function Documentation

void SiTrivialInduceChargeOnStrips::induce ( const SiChargeCollectionDrifter::collection_type &  collection_points, const StripGeomDetUnit &  det, std::vector< float > &  localAmplitudes, size_t &  recordMinAffectedStrip, size_t &  recordMaxAffectedStrip, const TrackerTopology *  tTopo  ) const

virtual

Implements SiInduceChargeOnStrips.

Definition at line 90 of file SiTrivialInduceChargeOnStrips.cc.

References induceVector().

                                           {


   induceVector(collection_points, det, localAmplitudes, recordMinAffectedStrip, recordMaxAffectedStrip, tTopo);

  /*
  auto ominA=recordMinAffectedStrip, omaxA=recordMaxAffectedStrip;
  std::vector<float> oampl(localAmplitudes);    
  induceOriginal(collection_points, det, oampl, ominA, omaxA, tTopo);

  //  std::cout << "orig " << ominA << " " << omaxA << " ";
  //for (auto a : oampl) std::cout << a << ",";
  //std::cout << std::endl;

  auto minA=recordMinAffectedStrip, maxA=recordMaxAffectedStrip;
  std::vector<float> ampl(localAmplitudes);
  induceVector(collection_points, det, ampl, minA, maxA, tTopo);

  // std::cout << "vect " << minA << " " << maxA << " ";          
  //for (auto a :   ampl) std::cout << a << ",";
  //std::cout << std::endl;
 
  float diff=0;
  for (size_t i=0; i!=ampl.size(); ++i) { diff = std::max(diff,ampl[i]>0 ? std::abs(ampl[i]-oampl[i])/ampl[i] : 0);}
  if (diff> 1.e-4) {
    std::cout << diff << std::endl;
    std::cout << "orig " << ominA << " " << omaxA << " ";
//    for (auto a : oampl) std::cout << a << ",";
    std::cout << std::endl;
    std::cout << "vect " << minA << " " << maxA << " ";
//    for (auto a : ampl) std::cout << a << ",";
    std::cout << std::endl;
  }

  localAmplitudes.swap(ampl);
  recordMinAffectedStrip=minA;
  recordMaxAffectedStrip=maxA;
  */

}

void SiTrivialInduceChargeOnStrips::induceOriginal ( const SiChargeCollectionDrifter::collection_type &  collection_points, const StripGeomDetUnit &  det, std::vector< float > &  localAmplitudes, size_t &  recordMinAffectedStrip, size_t &  recordMaxAffectedStrip, const TrackerTopology *  tTopo  ) const

private

Definition at line 250 of file SiTrivialInduceChargeOnStrips.cc.

References funct::abs(), prof2calltree::count, geVperElectron, StripTopology::localPitch(), bookConverter::max, min(), Nsigma, StripTopology::nstrips(), StripGeomDetUnit::specificTopology(), StripTopology::strip(), and ecaldqm::topology().

                                               {


  auto const & coupling = signalCoupling[typeOf(det,tTopo)];
  const StripTopology& topology = dynamic_cast<const StripTopology&>(det.specificTopology());
  size_t Nstrips =  topology.nstrips();

  if (!collection_points.empty()) count.dep(collection_points.size());

  for (auto signalpoint = collection_points.begin();  signalpoint != collection_points.end();  signalpoint++ ) {
    
    //In strip coordinates:
    double chargePosition = topology.strip(signalpoint->position());
    double chargeSpread = signalpoint->sigma() / topology.localPitch(signalpoint->position());
    
    size_t fromStrip  = size_t(std::max( 0,          int(std::floor( chargePosition - Nsigma*chargeSpread))));
    size_t untilStrip = size_t(std::min( Nstrips, size_t(std::ceil( chargePosition + Nsigma*chargeSpread) )));

    count.str(std::max(0,int(untilStrip)-int(fromStrip)));
    for (size_t strip = fromStrip;  strip < untilStrip; strip++) {

      double chargeDepositedOnStrip = chargeDeposited( strip, Nstrips, signalpoint->amplitude() / geVperElectron, chargeSpread, chargePosition);

      size_t affectedFromStrip  = size_t(std::max( 0, int(strip - coupling.size() + 1)));
      size_t affectedUntilStrip = size_t(std::min( Nstrips, strip + coupling.size())   );  
      for (size_t affectedStrip = affectedFromStrip;  affectedStrip < affectedUntilStrip;  affectedStrip++) {
    localAmplitudes.at( affectedStrip ) += chargeDepositedOnStrip * coupling.at(abs( affectedStrip - strip )) ;
      }

      if( affectedFromStrip  < recordMinAffectedStrip ) recordMinAffectedStrip = affectedFromStrip;
      if( affectedUntilStrip > recordMaxAffectedStrip ) recordMaxAffectedStrip = affectedUntilStrip;
    }
  }

}

void SiTrivialInduceChargeOnStrips::induceVector ( const SiChargeCollectionDrifter::collection_type &  collection_points, const StripGeomDetUnit &  det, std::vector< float > &  localAmplitudes, size_t &  recordMinAffectedStrip, size_t &  recordMaxAffectedStrip, const TrackerTopology *  tTopo  ) const

private

include last strip

do crosstalk... (can be done better, most probably not worth)

Definition at line 138 of file SiTrivialInduceChargeOnStrips.cc.

References funct::abs(), approx_erf(), DeDxDiscriminatorTools::charge(), constexpr, prof2calltree::count, delta, f, geVperElectron, i, j, gen::k, GetRecoTauVFromDQM_MC_cff::kk, StripTopology::localPitch(), bookConverter::max, min(), N, Nsigma, StripTopology::nstrips(), position, StripGeomDetUnit::specificTopology(), mathSSE::sqrt(), StripTopology::strip(), ecaldqm::topology(), and relativeConstraints::value.

Referenced by induce().

                                               {


  auto const & coupling = signalCoupling[typeOf(det,tTopo)];
  const StripTopology& topology = dynamic_cast<const StripTopology&>(det.specificTopology());
  const int Nstrips =  topology.nstrips();

  if (Nstrips == 0) return;

  const int NP = collection_points.size();
  if(0==NP) return;

  constexpr int MaxN = 512;
  // if NP too large split...

  for (int ip=0; ip<NP; ip+=MaxN) {
    auto N = std::min(NP-ip,MaxN);

    count.dep(N);
    float amplitude[N];
    float chargePosition[N];
    float chargeSpread[N];
    int fromStrip[N];
    int nStrip[N];

    // load not vectorize
    //In strip coordinates:
    for (int i=0; i!=N;++i) {
      auto j = ip+i;
      if (0==collection_points[j].amplitude()) count.zero();
      chargePosition[i]=topology.strip(collection_points[j].position());
      chargeSpread[i]= collection_points[j].sigma() / topology.localPitch(collection_points[j].position());
      amplitude[i]=0.5f*collection_points[j].amplitude() / geVperElectron;
    }
    
    // this vectorize
    for (int i=0; i!=N;++i) {
      fromStrip[i]  = std::max( 0,  int(std::floor( chargePosition[i] - Nsigma*chargeSpread[i])) );
      nStrip[i] = std::min( Nstrips, int(std::ceil( chargePosition[i] + Nsigma*chargeSpread[i])) ) - fromStrip[i];
    }
    int tot=0;
    for (int i=0; i!=N;++i) tot += nStrip[i];
    tot+=N; // add last strip 
    count.val(tot);
    float value[tot];
    
    // assign relative position (lower bound of strip) in value;
    int kk=0;
    for (int i=0; i!=N;++i) {
      auto delta = 1.f/(std::sqrt(2.f)*chargeSpread[i]);
      auto pos = delta*(float(fromStrip[i])-chargePosition[i]);
      for (int j=0;j<=nStrip[i]; ++j)  
    value[kk++] = pos+float(j)*delta;  
    }
    assert(kk==tot);
    
    // main loop fully vectorized
    for (int k=0;k!=tot; ++k)
      value[k] = approx_erf(value[k]);
    
    // saturate 0 & NStrips strip to 0 and 1???
    kk=0;
    for (int i=0; i!=N;++i) {
      if (0 == fromStrip[i])  value[kk]=0;
      kk+=nStrip[i];
      if (Nstrips == fromStrip[i]+nStrip[i]) value[kk]=1.f;
      ++kk;
    }
    assert(kk==tot);
    
    // compute integral over strip (lower bound becomes the value)
    for (int k=0;k!=tot-1; ++k)
      value[k]-=value[k+1];  // this is negative!
    
    
    float charge[Nstrips]; for (int i=0;i!=Nstrips; ++i) charge[i]=0;
    kk=0;
    for (int i=0; i!=N;++i){ 
      for (int j=0;j!=nStrip[i]; ++j)
    charge[fromStrip[i]+j]-= amplitude[i]*value[kk++];
      ++kk; // skip last "strip"
    }
    assert(kk==tot);
    
    
    int minA=recordMinAffectedStrip, maxA=recordMaxAffectedStrip;
    int sc = coupling.size();
    for (int i=0;i!=Nstrips; ++i) {
      int strip = i;
      if (0==charge[i]) continue;
      auto affectedFromStrip  = std::max( 0, strip - sc + 1);
      auto affectedUntilStrip = std::min(Nstrips, strip + sc);  
      for (auto affectedStrip=affectedFromStrip;  affectedStrip < affectedUntilStrip;  ++affectedStrip)
    localAmplitudes[affectedStrip] += charge[i] * coupling[std::abs(affectedStrip - strip)] ;
      
      if( affectedFromStrip  < minA ) minA = affectedFromStrip;
      if( affectedUntilStrip > maxA ) maxA = affectedUntilStrip;
    }
    recordMinAffectedStrip=minA;
    recordMaxAffectedStrip=maxA;
  }  // end loop ip

}

Member Data Documentation

const float SiTrivialInduceChargeOnStrips::geVperElectron

private

Definition at line 41 of file SiTrivialInduceChargeOnStrips.h.

Referenced by induceOriginal(), and induceVector().

const float SiTrivialInduceChargeOnStrips::Nsigma

private

Definition at line 40 of file SiTrivialInduceChargeOnStrips.h.

Referenced by induceOriginal(), and induceVector().

const std::vector<std::vector<float> > SiTrivialInduceChargeOnStrips::signalCoupling

private

Definition at line 38 of file SiTrivialInduceChargeOnStrips.h.