#include <MaterialBudgetCastorHistos.h>

Public Member Functions
void	fillEndTrack ()
void	fillPerStep (const G4Step *)
void	fillStartTrack (const G4Track *)
	MaterialBudgetCastorHistos (const edm::ParameterSet &p)
virtual	~MaterialBudgetCastorHistos ()
Private Member Functions
void	book ()
void	fillHisto (int id, int ix)
Private Attributes
int	binEta
int	binPhi
double	eta
double	etaHigh
double	etaLow
bool	fillHistos
int	id1
int	id2
double	intLen
std::vector< double >	intLength
std::vector< std::string >	matList
TProfile *	me100 [maxSet]
TProfile2D *	me1000 [maxSet]
TProfile2D *	me1100 [maxSet]
TH2F *	me1200 [maxSet]
TProfile *	me200 [maxSet]
TProfile *	me300 [maxSet]
TH1F *	me400 [maxSet]
TProfile *	me500 [maxSet]
TProfile *	me600 [maxSet]
TProfile *	me700 [maxSet]
TH1F *	me800 [maxSet]
TProfile2D *	me900 [maxSet]
double	phi
bool	printSum
double	radLen
std::vector< double >	radLength
double	stepLen
std::vector< double >	stepLength
int	steps
Static Private Attributes
static const int	maxSet = 20

Detailed Description

Definition at line 17 of file MaterialBudgetCastorHistos.h.

Constructor & Destructor Documentation

MaterialBudgetCastorHistos::MaterialBudgetCastorHistos ( const edm::ParameterSet & p )

Definition at line 15 of file MaterialBudgetCastorHistos.cc.

References binEta, binPhi, book(), etaHigh, etaLow, fillHistos, edm::ParameterSet::getUntrackedParameter(), pi, and printSum.

                                                                              {

  binEta      = p.getUntrackedParameter<int>("NBinEta", 100);
  binPhi      = p.getUntrackedParameter<int>("NBinPhi", 180);
  etaLow      = p.getUntrackedParameter<double>("EtaLow",  5.0);
  etaHigh     = p.getUntrackedParameter<double>("EtaHigh", 7.0);
  fillHistos  = p.getUntrackedParameter<bool>("FillHisto", true);
  printSum    = p.getUntrackedParameter<bool>("PrintSummary", false);
  edm::LogInfo("MaterialBudget") << "MaterialBudgetCastorHistos: FillHisto : "
                                 << fillHistos << " PrintSummary " << printSum
                                 << " == Eta plot: NX " << binEta << " Range "
                                 << etaLow << " (" << -etaHigh << ") : " 
                                 << etaHigh << " (" << -etaLow <<") Phi plot: "
                                 << "NX " << binPhi << " Range " << -pi << ":"
                                 << pi  << " (Eta limit " << etaLow << ":" 
                                 << etaHigh <<")";
  if (fillHistos) book();

}

MaterialBudgetCastorHistos::~MaterialBudgetCastorHistos ( ) [virtual]

Definition at line 35 of file MaterialBudgetCastorHistos.cc.

                                                        {
  edm::LogInfo("MaterialBudget") << "MaterialBudgetCastorHistos: Save user "
                                 << "histos ===";
}

Member Function Documentation

void MaterialBudgetCastorHistos::book ( ) [private]

Definition at line 169 of file MaterialBudgetCastorHistos.cc.

References binEta, binPhi, etaHigh, etaLow, i, edm::Service< T >::isAvailable(), maxEta, maxSet, me100, me1000, me1100, me1200, me200, me300, me400, me500, me600, me700, me800, me900, benchmark_cfg::minEta, mergeVDriftHistosByStation::name, pi, and indexGen::title.

Referenced by MaterialBudgetCastorHistos().

                                      {

  // Book histograms
  edm::Service<TFileService> tfile;
  
  if ( !tfile.isAvailable() )
    throw cms::Exception("BadConfig") << "TFileService unavailable: "
                                      << "please add it to config file";

  double maxPhi=pi;
  edm::LogInfo("MaterialBudget") << "MaterialBudgetCastorHistos: Booking user "
                                 << "histos === with " << binEta << " bins "
                                 << "in eta from " << etaLow << " to "
                                 << etaHigh << " and " << binPhi << " bins "
                                 << "in phi from " << -maxPhi << " to " 
                                 << maxPhi;
  
  char  name[10], title[80], tag1[10], tag2[15];
  // total X0
  for (int i=0; i<maxSet; i++) {
    double minEta=etaLow;
    double maxEta=etaHigh;
    int    ireg = i;
    if (i > 9) {
      minEta = -etaHigh;
      maxEta = -etaLow;
      ireg  -= 10;
      sprintf (tag2, " (-ve Eta Side)");
    } else {
      sprintf (tag2, " (+ve Eta Side)");
    }
    if ((i%2) == 0) {
      ireg  /= 2;
      sprintf (tag1, " == Start");
    } else {
      ireg   = (ireg-1)/2;
      sprintf (tag1, " == End");
    }
    sprintf(name, "%d", i+100);
    sprintf(title, "MB(X0) prof Eta in region %d%s%s", ireg, tag1, tag2);
    me100[i] =  tfile->make<TProfile>(name, title, binEta, minEta, maxEta);
    sprintf(name, "%d", i+200);
    sprintf(title, "MB(L0) prof Eta in region %d%s%s", ireg, tag1, tag2);
    me200[i] = tfile->make<TProfile>(name, title, binEta, minEta, maxEta);
    sprintf(name, "%d", i+300);
    sprintf(title, "MB(Step) prof Eta in region %d%s%s", ireg, tag1, tag2);
    me300[i] = tfile->make<TProfile>(name, title, binEta, minEta, maxEta);
    sprintf(name, "%d", i+400);
    sprintf(title, "Eta in region %d%s%s", ireg, tag1, tag2);
    me400[i] = tfile->make<TH1F>(name, title, binEta, minEta, maxEta);
    sprintf(name, "%d", i+500);
    sprintf(title, "MB(X0) prof Ph in region %d%s%s", ireg, tag1, tag2);
    me500[i] = tfile->make<TProfile>(name, title, binPhi, -maxPhi, maxPhi);
    sprintf(name, "%d", i+600);
    sprintf(title, "MB(L0) prof Ph in region %d%s%s", ireg, tag1, tag2);
    me600[i] = tfile->make<TProfile>(name, title, binPhi, -maxPhi, maxPhi);
    sprintf(name, "%d", i+700);
    sprintf(title, "MB(Step) prof Ph in region %d%s%s", ireg, tag1, tag2);
    me700[i] = tfile->make<TProfile>(name, title, binPhi, -maxPhi, maxPhi);
    sprintf(name, "%d", i+800);
    sprintf(title, "Phi in region %d%s%s", ireg, tag1, tag2);
    me800[i] = tfile->make<TH1F>(name, title, binPhi, -maxPhi, maxPhi);
    sprintf(name, "%d", i+900);
    sprintf(title, "MB(X0) prof Eta Phi in region %d%s%s", ireg, tag1, tag2);
    me900[i] = tfile->make<TProfile2D>(name, title, binEta/2, minEta, maxEta,
                                       binPhi/2, -maxPhi, maxPhi);
    sprintf(name, "%d", i+1000);
    sprintf(title, "MB(L0) prof Eta Phi in region %d%s%s", ireg, tag1, tag2);
    me1000[i]= tfile->make<TProfile2D>(name, title, binEta/2, minEta, maxEta,
                                       binPhi/2, -maxPhi, maxPhi);
    sprintf(name, "%d", i+1100);
    sprintf(title, "MB(Step) prof Eta Phi in region %d%s%s", ireg, tag1, tag2);
    me1100[i]= tfile->make<TProfile2D>(name, title, binEta/2, minEta, maxEta,
                                       binPhi/2, -maxPhi, maxPhi);
    sprintf(name, "%d", i+1200);
    sprintf(title, "Eta vs Phi in region %d%s%s", ireg, tag1, tag2);
    me1200[i]= tfile->make<TH2F>(name, title, binEta/2, minEta, maxEta, 
                                 binPhi/2, -maxPhi, maxPhi);
  }

  edm::LogInfo("MaterialBudget") << "MaterialBudgetCastorHistos: Booking user "
                                 << "histos done ===";

}

void MaterialBudgetCastorHistos::fillEndTrack ( )

Definition at line 154 of file MaterialBudgetCastorHistos.cc.

References fillHisto(), fillHistos, id1, id2, intLength, matList, printSum, radLength, and stepLength.

Referenced by MaterialBudgetHcal::update().

                                              {

  if (fillHistos) {
    if (id1 != 0) fillHisto(id1, 1);
    if (id2 != 0) fillHisto(id2, 1);
  }
  if (printSum) {
    for (unsigned int ii=0; ii<matList.size(); ii++) {
      edm::LogInfo("MaterialBudget") << matList[ii] << "\t" << stepLength[ii]
                                     << "\t" << radLength[ii] << "\t"
                                     << intLength[ii];
    }
  }
}

void MaterialBudgetCastorHistos::fillHisto	(	int	id,
		int	ix
	)		`[private]`

Definition at line 254 of file MaterialBudgetCastorHistos.cc.

References eta, etaHigh, etaLow, intLen, LogDebug, me100, me1000, me1100, me1200, me200, me300, me400, me500, me600, me700, me800, me900, phi, radLen, and stepLen.

Referenced by fillEndTrack(), and fillPerStep().

                                                         {

  int ii = 2*(id-1) + ix;
  double etaAbs = eta;
  if (eta < 0) {
    etaAbs = -eta;
    ii    += 10;
  }
  LogDebug("MaterialBudget") << "MaterialBudgetCastorHistos:FillHisto "
                             << "called with index " << id << ":" << ix 
                             << ":" << ii << " eta " << etaAbs << " (" 
                             << eta << ") integrated  step " << stepLen 
                             << " X0 " << radLen << " Lamda " << intLen;
  
  me100[ii]->Fill(eta, radLen);
  me200[ii]->Fill(eta, intLen);
  me300[ii]->Fill(eta, stepLen);
  me400[ii]->Fill(eta);

  if (etaAbs >= etaLow && etaAbs <= etaHigh) {
    me500[ii]->Fill(phi, radLen);
    me600[ii]->Fill(phi, intLen);
    me700[ii]->Fill(phi, stepLen);
    me800[ii]->Fill(phi);
  }

  me900[ii]->Fill(eta, phi, radLen);
  me1000[ii]->Fill(eta, phi, intLen);
  me1100[ii]->Fill(eta, phi, stepLen);
  me1200[ii]->Fill(eta, phi);
    
}

void MaterialBudgetCastorHistos::fillPerStep ( const G4Step * aStep )

Definition at line 70 of file MaterialBudgetCastorHistos.cc.

References fillHisto(), fillHistos, newFWLiteAna::found, g, id1, id2, intLen, intLength, testEve_cfg::level, LogDebug, matList, mergeVDriftHistosByStation::name, printSum, radLen, radLength, launcher::step, stepLen, and stepLength.

Referenced by MaterialBudgetHcal::update().

                                                                {

  G4Material * material = aStep->GetPreStepPoint()->GetMaterial();
  double step    = aStep->GetStepLength();
  double radl    = material->GetRadlen();
  double intl    = material->GetNuclearInterLength();
  double density = material->GetDensity() / (g/cm3);

  int    id1Old   = id1;
  int    id2Old   = id2;
  const G4VTouchable* touch = aStep->GetPreStepPoint()->GetTouchable();
  std::string         name  = touch->GetVolume(0)->GetName();
  std::string         matName = material->GetName();
  if (printSum) {
    bool found = false;
    for (unsigned int ii=0; ii<matList.size(); ii++) {
      if (matList[ii] == matName) {
        stepLength[ii] += step;
        radLength[ii]  += (step/radl);
        intLength[ii]  += (step/intl);
        found           = true;
        break;
      }
    }
    if (!found) {
      matList.push_back(matName);
      stepLength.push_back(step);
      radLength.push_back(step/radl);
      intLength.push_back(step/intl);
    }
    edm::LogInfo("MaterialBudget") << name << " " << step << " " << matName 
                                   << " " << stepLen << " " << step/radl << " " 
                                   << radLen << " " <<step/intl << " " <<intLen;
  } else {
    edm::LogInfo("MaterialBudget") << "MaterialBudgetCastorHistos: Step at " 
                                   << name << " Length " << step << " in " 
                                   << matName << " of density " << density 
                                   << " g/cc; Radiation Length " <<radl <<" mm;"
                                   << " Interaction Length " << intl << " mm\n"
                                   << "                          Position " 
                                   << aStep->GetPreStepPoint()->GetPosition()
                                   << " Cylindrical R "
                                   <<aStep->GetPreStepPoint()->GetPosition().perp()
                                   << " Length (so far) " << stepLen << " L/X0 "
                                   << step/radl << "/" << radLen << " L/Lambda "
                                   << step/intl << "/" << intLen;
  }

  int level = ((touch->GetHistoryDepth())+1);
  std::string name1="XXXX", name2="XXXX";
  if (level>3) name1 = touch->GetVolume(level-4)->GetName();
  if (level>4) name2 = touch->GetVolume(level-5)->GetName();
  if (name1 == "CAST") {
    id1 = 1;
    if      (name2 == "CAEC") id2 = 2;
    else if (name2 == "CAHC") id2 = 3;
    else if (name2 == "CEDC") id2 = 4;
    else if (name2 == "CHDC") id2 = 5;
    else                      id2 = 0;
  } else {
    id1 = id2 = 0;
  }
  LogDebug("MaterialBudget") << "MaterialBudgetCastorHistos: Level " << level
                             << " Volume " << name1 << " and " << name2
                             << " ID1 " << id1 << " (" << id1Old << ") ID2 "
                             << id2 << " (" << id2Old << ")";

  if (fillHistos) {
    if (id1 != id1Old) {
      if (id1 == 0) fillHisto(id1Old, 1);
      else          fillHisto(id1,    0);
    }
    if (id2 != id2Old) {
      if (id2 == 0) fillHisto(id2Old, 1);
      else          fillHisto(id2,    0);
    }
  }

  stepLen += step;
  radLen  += step/radl;
  intLen  += step/intl;
}

void MaterialBudgetCastorHistos::fillStartTrack ( const G4Track * aTrack )

Definition at line 40 of file MaterialBudgetCastorHistos.cc.

References dir, eta, id1, id2, intLen, intLength, matList, phi, printSum, radLen, radLength, stepLen, stepLength, and steps.

Referenced by MaterialBudgetHcal::update().

                                                                     {

  id1    = id2    = steps   = 0;
  radLen = intLen = stepLen = 0;

  const G4ThreeVector& dir = aTrack->GetMomentum() ;
  if (dir.theta() != 0 ) {
    eta = dir.eta();
  } else {
    eta = -99;
  }
  phi = dir.phi();
  double theEnergy = aTrack->GetTotalEnergy();
  int    theID     = (int)(aTrack->GetDefinition()->GetPDGEncoding());

  if (printSum) {
    matList.clear();
    stepLength.clear();
    radLength.clear();
    intLength.clear();
  }

  edm::LogInfo("MaterialBudget") << "MaterialBudgetCastorHistos: Track " 
                                 << aTrack->GetTrackID() << " Code " << theID
                                 << " Energy " << theEnergy/GeV << " GeV; Eta "
                                 << eta << " Phi " << phi/deg << " PT "
                                 << dir.perp()/GeV << " GeV *****";
}