DSTUtil::Calo Namespace Reference

Functions
int	GetHCALEnergy (const MIPPEventSummary *evt, float ptt, float hcE[2])
	Function returns the HCAL energy.
int	MakeTrkEMCalMatch (const MIPPEventSummary *evt, const MIPPTrackSummary *trk, float shXYEn[4], int shInd[2])
	Function makes track - emcal matching for the specified track If match is on then return the shower position and energy.
int	IsTrkHasCloseTrk (const MIPPEventSummary *evt, const MIPPTrackSummary *, int it)
	Function makes test for the presence of another track within 10 cm at EMCAL area.
int	GetEMCalHCalEnergies (float sumADC[2], float calE[2], int cid)
int	GetCalEnergy (const MIPPEventSummary *evt, float calE[3])
int	NumTrksHCAL (const MIPPEventSummary *evt)
	Function calculate the number of tracks in HCAL aperture (+/-55 cm).

---

Function Documentation

int DSTUtil::Calo::GetHCALEnergy	(	const MIPPEventSummary *	evt,
		float	ptt,
		float	hcE[2]
	)

Function returns the HCAL energy.

Definition at line 518 of file DSTUtil.cxx.

References MIPPEventSummary::GetHCal(), MIPPHCalSummary::hcalADC, RIMCConfig::kMCRunOffset, MIPPEventSummary::NHCal(), and MIPPEventSummary::run.

Referenced by NuMIPID::Fill(), fillCALOSummary(), and main().

{

  // HCAL ADC corrections made by using the minimum ionizing tracks (muons)
  // float corrADC[8]={0.77,1.43,0.54,1.28,1.53,1.72,0.78,0.87};
  // Apply the gain correction using 120 protons, where PMT's East - West 
  // balancing was made by cumulative approach
  float corrADC[8]={0.77515,1.2456,0.444,1.33006,1.1326,1.58409,1.2915,1.0058};
  float hcalsum=0.0; 
  // scale factor derived from 58 GeV/c and 120 GeV/c hadrons data 
  // Apply the HCAL energy scale dependable to the track momentum
  // energyscale = (16.5275 -0.0148231 * ptot)/1000.
  double enScale[2] = {16.5275, -0.0148231}; // for data and MC use
  double energyscale = 14.80357;   // default for 120 GeV/c protons

  energyscale = 14.80357;
  if(ptt<140.) energyscale = (enScale[0] + enScale[1]*ptt)/1000.;

 
  hcE[0] = 0.;
  hcE[1] = 0.;

  if(evt->NHCal()<1) return 9;    // there is no data

  MIPPHCalSummary* fHCal = evt->GetHCal(0);
  for(int j=0; j<8; ++j) {
    if(fHCal->hcalADC[j] > 0) {
      if (evt->run < RIMCConfig::kMCRunOffset) 
    hcalsum += (fHCal->hcalADC[j])*corrADC[j];
      else 
    hcalsum += fHCal->hcalADC[j];
    }
  } // end j loop
  hcE[0] = hcalsum * energyscale;    // GeV units
  hcE[1] = hcalsum;                  // raw ADC sum

  return 0;


} // end GetHCALEnergy function

int DSTUtil::Calo::MakeTrkEMCalMatch	(	const MIPPEventSummary *	evt,
		const MIPPTrackSummary *	trk,
		float	shXYEn[4],
		int	shInd[2]
	)

Function makes track - emcal matching for the specified track If match is on then return the shower position and energy.

Definition at line 564 of file DSTUtil.cxx.

References ds, ecalEventNum, MIPPEventSummary::ECalNumShwX, MIPPEventSummary::ECalNumShwY, MIPPTrackSummary::emcalx, MIPPEventSummary::evtnum, MIPPEventSummary::GetCalo(), MIPPEventSummary::NCalo(), numshx, numshy, MIPPTrackSummary::ptot, shenx, sheny, shx, and shy.

Referenced by NuMIPID::Fill(), fillCALOSummary(), main(), and RICHPid::MakeWeightedHists().

{

  // scale factor derived from 58 GeV/c and 120 GeV/c hadron data 
  double enEmcal[2] = {0.33083,  0.001748061}; // valid for data and MC
  double emscfactor=0.53413; // default for 120 GeV/c protons (in MeV)  
  emscfactor = (enEmcal[0] + enEmcal[1]*trk->ptot)/1000.;  // now in GeV
  //float emcalSum=0.0;
  float pmatchWin[6] = {8., 6., 4., 3., 2., 1.};
  float dxx,dyy,dss[2]={5.0, 6.0};

  if(evt->evtnum != ecalEventNum) {  // is it same event?
    ecalEventNum = evt->evtnum;
// use codes below only for the first call in event
    if(evt->NCalo() < 1) {
      numshx = 0;
      numshy = 0;
      shXYEn[0] = 200.;
      shXYEn[1] = 210.;
      shXYEn[2] = 0.;
      shXYEn[3] = 0.;
      shInd[0]  = -1;
      shInd[1]  = -2;
      return 9;
    } // end if(evt->NCalo() < 1) statement

    numshx = evt->ECalNumShwX;
    numshy = evt->ECalNumShwY;
    for(int ical=0; ical < evt->NCalo(); ++ical) {
      MIPPCaloSummary* fcal = evt->GetCalo(ical);
      shx[ical] = fcal->EShXP;
      shy[ical] = fcal->EShYP;
      shenx[ical] = fcal->EShXEn;
      sheny[ical] = fcal->EShYEn;
    } // end ical loop
  } // end if(evt->evtnum != ecalEventNum) statement

  if(numshx==0 || numshy==0) {  // do not make match if one of them is 0
    shXYEn[0] = 200.;
    shXYEn[1] = 210.;
    shXYEn[2] = 0.;
    shXYEn[3] = 0.;
    shInd[0]  = -1;
    shInd[1]  = -2;
    return 9;
  } // end if(numshx==0 || numshy==0) statement

  // now start actual matching
  int jmom,jx,jy;
  float ds, dx, dy;
  jmom = int(trk->ptot/10.0);
  if(jmom < 0) jmom = 0;
  if(jmom > 5) jmom = 5;
  jx = -1;
  jy = -1;
  // apply momentum dependence for the matching window
  dxx = dss[0] + pmatchWin[jmom];
  ds = 100.;
  for(int ix=0; ix<numshx; ++ix) {
    dx = fabs(trk->emcalx[0]-shx[ix]);
    if(dx<ds && dx<dxx) { // match is on
      jx = ix;
      ds = dx;
    } // end if(dx<ds && dx<dxx) test
  } // end ix loop
  if(ds<dxx && jx>-1) { // closest match in X is found

    // now in Y-view
    dyy = dss[1] + pmatchWin[jmom];
    ds = 100.;
    for(int iy=0; iy<numshy; ++iy) {
      dy = fabs(trk->emcalx[1]-shy[iy]);
      if(dy<dyy && dy<ds) {  //match is on
    jy = iy;
    ds = dy;
      } // end match test
    } // end iy loop
    if(ds<dyy && jy>-1) { // closest match in Y view is found

      // now place positions and calculate an energy in GeV
      shXYEn[0] = shx[jx];
      shXYEn[1] = shy[jy];
      float shC = shenx[jx] + sheny[jy];
      shXYEn[2] = emscfactor*shC;
      shXYEn[3] = shC;      // raw ADC sum
      shInd[0]  = jx;
      shInd[1]  = jy;

    } // end if(ds<dyy && jy>-1) statement
  } // end if(ds<dxx && jx>-1) statement

  if(jx == -1 || jy == -1) {
    shXYEn[0] = 200.;
    shXYEn[1] = 210.;
    shXYEn[2] = 0.;
    shXYEn[3] = 0.;
    shInd[0]  = -1;
    shInd[1]  = -2;
    return 9;
  } // if(jx == -1 || jy == -1) statement

  return 0;  // match is found, variables are filled above

} // end MakeTrkEMCalMatch function

int DSTUtil::Calo::IsTrkHasCloseTrk	(	const MIPPEventSummary *	evt,
		const MIPPTrackSummary *	trk,
		int	it
	)

Function makes test for the presence of another track within 10 cm at EMCAL area.

Definition at line 674 of file DSTUtil.cxx.

References MIPPEventSummary::evtnum, MIPPEventSummary::GetTrk(), MIPPEventSummary::NTrk(), nTrks, nTrksAtEmcal, trkEventNum, xEmf, and yEmf.

Referenced by main().

{
 
  float howClose=10.0, xEmTrk[300], yEmTrk[300];
  if(evt->evtnum != trkEventNum) { 
    trkEventNum = evt->evtnum;
// use below codes only for the first call in event
    int nTrks = evt->NTrk();
    if(nTrks > 300) nTrks = 300; // what is a size of the track bank?
    if(nTrks < 2) return 0;      // nothing to do if nTrks < 2

    // Find how many tracks pointing to the EMCal area
    // apply loose cuts for tracks, do not apply RICH window cut
    nTrksAtEmcal = 0;
    int wc4, wc5, wc6;
    for(int i=0; i<nTrks; ++i) {
      MIPPTrackSummary* trp = evt->GetTrk(i);
      xEmf[i] = -1;    // use -1 if track is not pointing to EMCal
      yEmf[i] = -1;
      xEmTrk[i] = 250.;
      yEmTrk[i] = 260.;

      if(trp->avgt < -30. || trp->avgt > 120.) continue;   // test for time
      if(trp->ntpchit < 20 || trp->ntpchit > 95) continue; // test space points
      wc4 = wc5 = wc6 = 0;
      if(trp->nwchits[3] > 0) wc4 = 1;
      if(trp->nwchits[4] > 0) wc5 = 1;
      if(trp->nwchits[5] > 0) wc6 = 1;
      if((wc4+wc5+wc6) < 1) continue;    // at least one downstream WC is on
      // position of track at EMCAL
      if(fabs(trp->emcalx[0])>85.0 || fabs(trp->emcalx[1])>85.0) continue;

      xEmf[i] = 0;         // track in EMCal area, reset flag to 0
      yEmf[i] = 0;
      xEmTrk[i] = trp->emcalx[0];
      yEmTrk[i] = trp->emcalx[1];
      nTrksAtEmcal++;             // count it
    } // end track loop

    if(nTrksAtEmcal < 2) {   // if less than 2 tracks, then nothing to compare
      return 0;
    } else {
      // now look in X-view
      for(int i1=0; i1<nTrks-1; ++i1) {  // outer loop
    if(xEmf[i1] == -1) continue;  // track not at EMCal
    for(int i2=i1+1; i2<nTrks; ++i2) { // inner loop
      if(xEmf[i2] == -1) continue;  // track not at EMCal
      if(fabs(xEmTrk[i1] - xEmTrk[i2])<howClose) { //there are close tracks
        xEmf[i1] = 5;  // set found flag for both tracks
        xEmf[i2] = 6;
      } // end if < 10 test
    } // end i2 inner loop
      } // end outer i1 loop

      // now in Y-view
      for(int i3=0; i3<nTrks-1; ++i3) {
    if(yEmf[i3] == -1) continue;
    for(int i4=i3+1; i4<nTrks; ++i4) {
      if(yEmf[i4] == -1) continue; 
      if(fabs(yEmTrk[i3] - yEmTrk[i4])<howClose) {
        yEmf[i3] = 7;  
        yEmf[i4] = 8;
      } // end if < 10 test
    } // end i4 loop
      } // end i3 loop

    } // end else statement

  } // end if(evt->evtnum != trkEventNum) test

  // usually come here for another calls in event 
  int iflag = 0;
  if(xEmf[it] > 1 || yEmf[it] > 1) iflag = 9;
  return iflag;  // iflag=9 if another close track found 
} // end IsTrkHasCloseTrk function

int DSTUtil::Calo::GetEMCalHCalEnergies	(	float	sumADC[2],
		float	calE[2],
		int	cid
	)

Definition at line 761 of file DSTUtil.cxx.

{
  float emcalSF[15]={0.37,0.37,0.37,0.37,0.37,0.37,0.37,0.438,0.489,0.0,
           0.41,0.47,0.0,0.43,0.52};
  float hcalSF[15]={10.57,10.57,10.57,10.57,10.57,10.57,10.57,15.05,15.01,0.0,
          15.38,15.26,0.0,15.68,14.95};

  calE[0] = 0.;
  calE[1] = 0.;
  if(cid<1 || (cid>3 && cid<8) || cid==13 || cid>15) return 9;
  if(sumADC[0]<0. || sumADC[1]<0.) return 8;
  calE[0] = sumADC[0] * emcalSF[cid-1]/1000.0;  // return in GeV
  calE[1] = sumADC[1] * hcalSF[cid-1]/1000.0;   // return in GeV
  return 0;
} // end GetEMCalHCalEnergies function

int DSTUtil::Calo::GetCalEnergy	(	const MIPPEventSummary *	evt,
		float	calE[3]
	)

Definition at line 778 of file DSTUtil.cxx.

References MIPPEventSummary::ECalNumShwX, MIPPEventSummary::ECalNumShwY, MIPPCaloSummary::EShXEn, MIPPEventSummary::GetCalo(), MIPPEventSummary::GetHCal(), MIPPHCalSummary::hcalADC, RIMCConfig::kMCRunOffset, MIPPEventSummary::NCalo(), MIPPEventSummary::NHCal(), nShx, nShy, and MIPPEventSummary::run.

Referenced by RICHPid::MakeWeightedHists(), and ValiCham::OneTrackEvent().

{

// get H & EM Cal Energies based on the Cal Summaries.
// energies are stored in the calE array
// calE[0]  = total HCal
// calE[1] = total EMCal
// calE[2] = total EM+H

  calE[0] = 0.;
  calE[1] = 0.;
  calE[2] = 0.;

  // don't correct hcal adc's if MC 
  bool mcRun = false; 
  if (evt->run >= RIMCConfig::kMCRunOffset) mcRun = true;

  // get hcal energy
  float sumhcal=0.0, sumHcal=0.0, sfh=14.80357, sfe=0.53413;
  float corrADC[8]={0.77515,1.2456,0.444,1.33006,1.1326,1.58409,1.2915,1.0058};

  for(int ical=0; ical<evt->NHCal(); ++ical) {
    if (ical > 0) continue; 
    MIPPHCalSummary* fHCal = evt->GetHCal(ical);
    for(int j=0; j<8; ++j) {
      if (mcRun) 
        sumhcal += fHCal->hcalADC[j];
      else 
        sumhcal += (fHCal->hcalADC[j]*corrADC[j]);
    } // end j loop
  } // end ical loop
  sumHcal = sumhcal/1000.;
  calE[0] = sumHcal*sfh;

  // get emcal energy 
  int nShx = evt->ECalNumShwX;
  int nShy = evt->ECalNumShwY;
  float emXen = 0.0, emYen = 0.0;
  
  for(int ical=0; ical < evt->NCalo(); ++ical) {
    MIPPCaloSummary* fcal = evt->GetCalo(ical);
    if(ical < nShx) emXen += fcal->EShXEn;
    if(ical < nShy) emYen += fcal->EShYEn;
  } // end ical loop
  calE[1] = (emXen + emYen)*sfe/1000.; // total EMCal Energy in GeV
  if (calE[1] < 0.) calE[1] = 0.0;  
 
  // total cal energy 
  calE[2] = calE[0] + calE[1]; // total CAL energy in GeV
  
  return 0;
}

int DSTUtil::Calo::NumTrksHCAL

(

const MIPPEventSummary *

evt

)

Function calculate the number of tracks in HCAL aperture (+/-55 cm).

Definition at line 833 of file DSTUtil.cxx.

References MIPPEventSummary::GetTrk(), MIPPEventSummary::NTrk(), and nTrks.

Referenced by fillCALOSummary(), and main().

{
  int nTrks = 0, wc4,wc5,wc6;
  for (int i = 0; i < evt->NTrk(); ++i) {
    MIPPTrackSummary* trk = evt->GetTrk(i);
    if(trk->avgt < -30. || trk->avgt > 120.) continue;   // test for time
    if(trk->ntpchit < 20 || trk->ntpchit > 95) continue; // test space points
    wc4 = wc5 = wc6 = 0;
    if(trk->nwchits[3] > 0) wc4 = 1;
    if(trk->nwchits[4] > 0) wc5 = 1;
    if(trk->nwchits[5] > 0) wc6 = 1;
    if((wc4+wc5+wc6) < 1) continue;    // at least one downstream WC is on
    // position of track at HCAL
    if(fabs(trk->hcalx[0])>55.0 || fabs(trk->hcalx[1])>55.0) continue;
    nTrks++;
  }
  return nTrks;

} // end NumTrksHCAL function

---