UtilStrip.cxx

Go to the documentation of this file.

#include "TMath.h"

#include "MessageService/MsgService.h"

#include "AtNuEvent/AtmosEvent.h"
#include "AtNuEvent/AtmosScintHit.h"
#include "AtNuEvent/AtmosStrip.h"
#include "AtNuEvent/AtmosMC.h"

#include "AtNuUtils/UtilStrip.h"

CVSID("$Id: UtilStrip.cxx,v 1.2 2006/07/19 22:41:29 gmieg Exp $");

const double C45 = TMath::Sqrt(2.) / 2.; //Cosine of 45 degrees

std::vector<int> UtilStrip::GangedStripsNegative(int strip) {
  std::vector<int> Gang;
  if(strip<0 || strip>191) return Gang;
  Gang.push_back(strip);

  int AddStrip = strip;
  //Climb Up
  while (true) {
    if(AddStrip%24 == 0) AddStrip += 47;
    else AddStrip += 23;
    if(AddStrip > 191) break;
    Gang.push_back(AddStrip);
  }

  //Climb Down
  AddStrip = strip;
  while (true) {
    if(AddStrip<47) AddStrip -= 23;
    else if((AddStrip-47)%24 == 0) AddStrip -= 47;
    else AddStrip -= 23;
    if(AddStrip < 0) break;
    Gang.push_back(AddStrip);
  }

  return Gang;
}

std::vector<int> UtilStrip::GangedStripsPositive(int strip) {
  std::vector<int> Gang;
  if(strip<0 || strip>191) return Gang;
  Gang.push_back(strip);

  int AddStrip = strip;
  while (true) {
    AddStrip += 24;
    if(AddStrip > 191) AddStrip -= 192;
    if(AddStrip == strip) break;
    Gang.push_back(AddStrip);
  }

  return Gang;
}

double UtilStrip::SumEnergy(std::set<SetStrip> strips) {
  double sumenergy = 0.0;
  std::set<SetStrip>::iterator itr = strips.begin();
  for(; itr!=strips.end(); itr++) sumenergy += itr->Energy;
  return sumenergy;
}

int UtilStrip::NStrips(std::set<SetStrip> strips) {
  int nstrips = 0;
  std::set<SetStrip>::iterator itr = strips.begin();
  for(; itr!=strips.end(); itr++) if(itr->Energy>0) nstrips++;
  return nstrips;
}

int UtilStrip::NDigits(std::set<SetStrip> strips) {
  int ndigits = 0;
  std::set<SetStrip>::iterator itr = strips.begin();
  for(; itr!=strips.end(); itr++)
    if(itr->Energy>0 && itr->DoubleEnded) ndigits++;
  return ndigits;
}

/*
 * SetStrip class
 */
SetStrip::SetStrip(int i, const AtmosStrip *strip) {
  TCIndex = i;

#ifdef PECORR
  Energy = strip->QPEcorr[0] + strip->QPEcorr[1];
#else
  Energy = strip->QPE[0] + strip->QPE[1];
#endif

  DoubleEnded = (bool)(strip->Ndigits - 1);
}

/*
 * StripCluster class
 */
void StripCluster::Clear() {
  NStrips = 0;
  SumU = SumV = SumZ = SumSqU = SumSqV = SumSqZ = 0.0;

  SumPE = 0.0;
  WSumU = WSumV = WSumZ = WSumSqU = WSumSqV = WSumSqZ = 0.0;

  Strips.clear();
}

void StripCluster::AddStripDList(const StripDList &sdl) {
  NStrips++;
  SumPE += sdl.PE;

  SumU += sdl.U;
  SumV += sdl.V;
  SumZ += sdl.Z;
  
  SumSqU += sdl.U * sdl.U;
  SumSqV += sdl.V * sdl.V;
  SumSqZ += sdl.Z * sdl.Z;

  WSumU += sdl.U * sdl.PE;
  WSumV += sdl.V * sdl.PE;
  WSumZ += sdl.Z * sdl.PE;
  
  WSumSqU += sdl.U * sdl.U * sdl.PE;
  WSumSqV += sdl.V * sdl.V * sdl.PE;
  WSumSqZ += sdl.Z * sdl.Z * sdl.PE;

  Strips.insert(SetStrip(sdl.TCIndex,false,sdl.PE));
}

double StripCluster::RMSU() {
  return TMath::Sqrt(SumSqU/(double)NStrips + (CentU()*CentU()));
}

double StripCluster::RMSV() {
  return TMath::Sqrt(SumSqV/(double)NStrips + (CentV()*CentV()));
}

double StripCluster::RMSZ() {
  return TMath::Sqrt(SumSqZ/(double)NStrips + (CentZ()*CentZ()));
}

StripDList::StripDList(int i) {
  TCIndex = i;
  InCluster = false;
  ClusterCount = 0;
  U = V = Z = 0.0;
}

bool StripDList::SetInCluster() {
  ClusterCount++;
  if(InCluster) return(true);
  InCluster = true;
  return(false);
}

bool operator<(const StripDList &sdl1, const StripDList &sdl2) {
  if(sdl1.VecDStrip.size() == 0 ||
     sdl2.VecDStrip.size() == 0) return false;
  return(sdl1.VecDStrip[0].Disp < sdl2.VecDStrip[0].Disp);
}

bool operator==(const StripDList &sdl1, const StripDList &sdl2) {
  return(sdl1.TCIndex == sdl2.TCIndex);
}

bool operator!=(const StripDList &sdl1, const StripDList &sdl2) {
  return(sdl1.TCIndex != sdl2.TCIndex);
}

---