|
//____________________________________________________________________ // // // //____________________________________________________________________ // // $Id: BrMultRdoModule.cxx,v 1.20 2002/05/08 16:30:55 cholm Exp $ // $Author: cholm $ // $Date: 2002/05/08 16:30:55 $ // $Copyright: (C) 2001 BRAHMS Collaboration <brahmlib@rhic.bnl.gov> // event #ifndef BRAT_BrMultRdoModule #include "BrMultRdoModule.h" #endif #ifndef ROOT_TDirectory #include "TDirectory.h" #endif #ifndef BRAT_BrRunInfoManager #include "BrRunInfoManager.h" #endif #ifndef WIN32 #include <iostream> #include <iomanip> #include <climits> #include <cfloat> #else #include <iostream.h> #include <iomanip.h> #include <limits.h> #include <float.h> #endif //____________________________________________________________________ ClassImp(BrMultRdoModule); //____________________________________________________________________ BrMultRdoModule::BrMultRdoModule() { // Default constructor. DO NOT USE SetState(kSetup); SetThresholdFactor(); SetSaturationChannel(); //Set OutlierMethod to kHitCorrection for Fraction of hit centrality // for kNoCorrection for Kansas multiplicity SetOutlierMethod(); SetAdcGapLimit(); fCalibration = 0; fParameters = 0; fEtaFills = 0; } //____________________________________________________________________ BrMultRdoModule::BrMultRdoModule(const Char_t* name, const Char_t* title) : BrModule(name, title), BrMultUtil(name) { // Named Constructor SetState(kSetup); SetThresholdFactor(); SetSaturationChannel(); //Set OutlierMethod to kHitCorrection for Fraction of hit centrality // for kNoCorrection for Kansas multiplicity SetOutlierMethod(); SetAdcGapLimit(); fCalibration = 0; fParameters = 0; fEtaFills = 0; } //____________________________________________________________________ void BrMultRdoModule::DefineHistograms() { // Define some histograms: // Dimension Contents // -------------------------------------------------- // 2D ADC vs Single // 2D Energy vs Single // 2D <Energy> vs Ring // 1D Total <Energy> // 2D Hits vs Ring // 1D Hits if (GetState() != kInit) { Stop("DefineHistograms", "Must be called after Init"); return; } if (!fParameters) { Failure("DefineHistograms", "no detector parameters"); return; } // Get some parameters Int_t nRows = fParameters->GetNoRows(); Int_t nRings = fParameters->GetNoRings(); Int_t nSingles = fParameters->GetNoModules(); TDirectory* saveDir = gDirectory; TDirectory* histDir = gDirectory->mkdir(GetName()); histDir->cd(); // Make histograms here fSingleAdcHisto = new TH2F("singleAdc", "Single vs ADC", nSingles+1, .5, nSingles + 1.5, (Int_t) ((fSingleAdcHigh-fSingleAdcLow ) / fSingleAdcComp), fSingleAdcLow, fSingleAdcHigh); fSingleCalAdcHisto = new TH2F("singleCalAdc", "Single vs Calibrated ADC", nSingles+1, .5, nSingles + 1.5, 421, -20, fBaseCalAdc); fSingleEnergyHisto = new TH2F("singleEnergy", "Single vs deposited Energy", nSingles+1, .5, nSingles + 1.5, 100, -0.05*fBaseEnergy, fBaseEnergy); fSingleMultHisto = new TH2F("singleMult", "Single vs Multiplicity", nSingles+1, .5, nSingles + 1.5, 100, 0, fBaseMult); fRingCalAdcHisto = new TH2F("ringCalAdc", "Ring vs Calibrated ADC", nRings, .5, nRings + .5, 100, 0, nRows * fBaseCalAdc); fRingEnergyHisto = new TH2F("ringEnergy", "Ring vs deposited energy", nRings, .5, nRings + .5, 100, 0, nRows * fBaseEnergy); fRingMultHisto = new TH2F("ringMult", "Ring vs Multiplicity", nRings, .5, nRings + .5, 100, 0, nRows * fBaseMult); fRingHitsHisto = new TH2F("ringHits", "Ring vs Hits", nRings, .5, nRings + .5, nRows, -.5, nRows + .5); fArrayCalAdcHisto = new TH1F("arrayCalAdc", "Array Calibrated ADC", 100, 0, nSingles * fBaseCalAdc); fArrayEnergyHisto = new TH1F("arrayEnergy", "Array deposited energy", 100, 0, nSingles * fBaseEnergy); fArrayMultHisto = new TH1F("arrayMult", "Array Multiplicity", 100, 0, nSingles * fBaseMult); fArrayHitsHisto = new TH1F("arrayHits", "Array Hits", nSingles+1, -.5, nSingles + 1.5); fVtxVsCalAdcHisto = new TH2F("vtxVsCalAdc", "V_{z} vs Calibrated ADC", 16, -100, 100, 100, 0, nSingles * fBaseCalAdc); fVtxVsMultHisto = new TH2F("vtzVsMult", "V_{z} vs Multiplicity", 16, -100, 100, 100, 0, nSingles * fBaseMult); fEtaDistHisto = new TH1F("etaDist", ""#frac{dN}{d#eta}"", 30, -3, 3); fdNdEtaSumHisto = new TH1F("dNdEtaSum", "#frac{dN}{d#eta} Sum", 30, -3, 3); fdNdEtaFillHisto = new TH1F("dNdEtaFill", "#frac{dN}{d#eta} normalisation", 30, -3, 3); fdNdEtaHisto = new TH1F("dNdEta", "#frac{dN}{d#eta}", 30, -3, 3); gDirectory = saveDir; } //____________________________________________________________________ void BrMultRdoModule::Begin() { fCurrentRun = BrRunInfoManager::Instance()->GetCurrentRun(); } //____________________________________________________________________ void BrMultRdoModule::Finish() { // Job-level finalisation // Normalise the dN/deta histogram to number of events and bin size SetState(kFinish); if (HistOn()) { fdNdEtaHisto->Divide(fdNdEtaSumHisto, fdNdEtaFillHisto); } } //____________________________________________________________________ void BrMultRdoModule::Print(Option_t* option) const { // Print module information // See BrModule::Print for options. // In addition this module defines the Option: // <fill in here> TString opt(option); opt.ToLower(); BrModule::Print(option); if (opt.Contains("d")) { cout << endl << " Original author: Christian Holm Christensen" << endl << " Last Modifications: " << endl << " $Author: cholm $" << endl << " $Date: 2002/05/08 16:30:55 $" << endl << " $Revision: 1.20 $ " << endl << endl << "-------------------------------------------------" << endl; cout << " Threshold Factor: " << fThresholdFactor << endl << " Saturation Chnl: " << fSaturationChannel << endl << " Outlier method: " << flush; switch(fOutlierMethod) { case kCloseNeighbors: cout << "Close Neighbors" << flush; break; case kFractionArray: cout << "Fraction of Array" << flush; break; default: cout << "No Correction" << flush; break; } cout << endl; } } //__________________________________________________________________ void BrMultRdoModule::CalibrateRings(BrMultRdo* rdo) { // PRIVATE METHOD: // Make the calibrated data for the individual detector groups i.e., // the rings. // This method does the outlier correction as well. There are three // ways of doing that: // // "Close Neighbors" // See the method BrMultRdoModule::CalibrateRingsCloseNeighbors() // // "Fraction Array" // See the method BrMultRdoModule::CalibrateRingsFractionArray() // // "None" // No correction for outliers, just sum it. // // As a side effect, this method also calculates the sum energy and // multiplicity over the full array, but a further correction need // to be applied. // Loop variables, must be decalered here so that MSVC++/KCC // doesn't get confused - poor sods, they should stick to // ANSI/ISO standard, alias they don't! Int_t i = 0; // Get some parameters on this detector Int_t nRings = fParameters->GetNoRings(); Int_t nSingles = fParameters->GetNoModules(); // First the outlier corrections if (fOutlierMethod == kCloseNeighbors) // "Close Neighbors" outlier correction CalibrateRingsCloseNeighbors(rdo); else if (fOutlierMethod == kFractionArray) // "Fraction Array" outlier correction CalibrateRingsFractionArray(rdo); else { // No outlier corrections // Zero the multitplicity in object for (i = 0; i < nRings; i++) { rdo->SetRingMultiplicity(i, 0); rdo->SetRingdNdEta(i, 0); rdo->SetRingAvgEta(i, 0); } float tot = 0; // Loop over the singles, and sum in each ring the multiplicty Int_t nHits[300]; for(i=0;i<300; i++) nHits[i] = 0; for (i = 0; i < nSingles; i++) { fCurrentRing = fCalibration->GetRingMap(i); fCurrentRingZ = fCalibration->GetRingPosition(fCurrentRing); if(rdo->GetSingleMultiplicity(i)>0) tot += rdo->GetSingleMultiplicity(i); // Check if single is on. if (fCurrentRing < 0 || fCurrentRing > nRings) // If not continue to next single continue; // Increment the ring if(fabs(fCurrentVtxZ)<30) ++nHits[fCurrentRing]; if(rdo->GetSingleMultiplicity(i)>-1&&rdo->GetSingleMultiplicity(i)<300.) { rdo->SetRingMultiplicity(fCurrentRing, rdo->GetRingMultiplicity(fCurrentRing) + rdo->GetSingleMultiplicity(i)); if (DebugLevel() > 30) cout << " (" << i << ") " << rdo->GetSingleMultiplicity(i) << " -> " << rdo->GetRingMultiplicity(fCurrentRing) << " ( " << fCurrentRing << " ) " << endl; if(fabs(fCurrentVtxZ)<30.) { Double_t ringScale; rdo->SetRingAvgEta(fCurrentRing, rdo->GetRingAvgEta(fCurrentRing) + CalibrateAvgEta(ringScale)); rdo->SetRingdNdEta(fCurrentRing, rdo->GetRingdNdEta(fCurrentRing) + ringScale* rdo->GetSingleMultiplicity(i)); } } } // Renormalize to get ring average if(fabs(fCurrentVtxZ)<30.) { for(i = 0; i< nRings; i++) { //Float_t hits = (Float_t) fCalibration->GetRingCnt(i); Float_t hits = (Float_t) nHits[i]; if(hits>0) { rdo->SetRingAvgEta(i, rdo->GetRingAvgEta(i)/hits); rdo->SetRingdNdEta(i, rdo->GetRingdNdEta(i)/hits); } else { rdo->SetRingAvgEta(i, -10.); rdo->SetRingdNdEta(i, 0.); } } } // Now correct for the number of singles hit. for (i = 0; i < nRings; i++) { if(fOutlierMethod == kHitCorrection) { if (rdo->GetRingHits(i) == 0) { rdo->SetRingMultiplicity(i, 0); rdo->SetRingEnergy(i, 0); continue; } rdo->SetRingMultiplicity(i, TMath::Max((Float_t) 0., rdo->GetRingMultiplicity(i) * nSingles / nRings / rdo->GetRingHits(i))); rdo->SetRingEnergy(i, TMath::Max((Float_t) 0., rdo->GetRingEnergy(i) * nSingles / nRings / rdo->GetRingHits(i))); } if (DebugLevel() > 30) cout << GetName() << " Ring " << i << "/" << nRings << " mult is " << rdo->GetRingMultiplicity(i) << " / " << tot << endl; // Increment the array sum rdo->SetArrayEnergy(rdo->GetArrayEnergy() + rdo->GetRingEnergy(i)); rdo->SetArrayMultiplicity(rdo->GetArrayMultiplicity() + rdo->GetRingMultiplicity(i)); } } // Fill histograms if (HistOn()) { for (i = 0; i < nRings; i++) { fRingEnergyHisto->Fill(i + 1, rdo->GetRingEnergy(i)); fRingMultHisto->Fill(i + 1, rdo->GetRingMultiplicity(i)); fEtaDistHisto->Fill(rdo->GetRingEta(i), rdo->GetRingMultiplicity(i)); fdNdEtaSumHisto->Fill(rdo->GetRingAvgEta(i), rdo->GetRingdNdEta(i)); fdNdEtaFillHisto->Fill(rdo->GetRingAvgEta(i)); } // Count the number of evenst that contribute to the dN/deta // histogram fEtaFills++; } } //__________________________________________________________________ void BrMultRdoModule::CalibrateRingsCloseNeighbors(BrMultRdo* rdo) { // PRIVATE METHOD: // Correct for outliers. // The "Close Neighbors" method calculates the avarage and variance // over a given ring and it's immediate neighbors. Then for each // single in the ring and it's immediate neighbors it tags a single as // an outlier, if it's energy signal is higher then the mean plus // three variances. For those singles tag as outliers, the method // substitute the energy with avarage over the the non-tagged singles // in the ring and it's immediate neighbors. // Loop variables, must be decalered here so that MSVC++/KCC // doesn't get confused - poor sods, they should stick to ANSI/ISO // standard, alias they don't! Int_t i = 0; Int_t j = 0; Int_t k = 0; // Get some parameters on this detector Int_t nRings = fParameters->GetNoRings(); Int_t nSingles = fParameters->GetNoModules(); for (i = 0; i < nRings; i++) { fCurrentRing = i; Int_t hits = 0; Double_t energyAvg = 0; Double_t energyVar = 0; Int_t accHits = 0; Double_t accEnergyAvg = 0; Double_t accMultAvg = 0; // Indicies of immediate neighbors Int_t lowRing = (i > 0 ? i - 1 : i); Int_t highRing = (i == nRings - 1 ? i : i + 1); // If the "Close Neighbors" method is used (fOutlierMethod == // kCloseNeighbors), we should calculate the avarage and variance // over the current ring, and it's immediate neighbors for (j = lowRing; j <= highRing; j++) { hits += rdo->GetRingHits(j); energyAvg += rdo->GetRingEnergy(j); // We used the member BrMultRdo::BrRing::fMultiplicity in the // method Br[Tile|Si]RdoModule::CalibrateIndividual to store the // sum square of the energy over the rings, in the member to // save some space. energyVar += rdo->GetRingMultiplicity(j); // cout << ringEnergy[j] << "/" << ringHits[j] << " " << flush; } // Find the average and variance over the current ring and it's // neighbors. energyAvg /= hits; energyVar = TMath::Sqrt(energyVar / hits - energyAvg * energyAvg); // cout << "Close neighbors <E>: " << energyAvg // << " s: " << energyVar // << " #: " << hits << endl; // For each single that contributed to any of the rings, check if // the energy is within the 3 variances from the mean. If so, // accept it, and calculate a new mean for (j = 0; j < nSingles; j++) { // Check if the single is part of the current ring or one of // it's close neighbors. That is, if the ring it belongs to is // in the interval [lowRing, highRing]. Note that we do a // negated test. if (lowRing > fCalibration->GetRingMap(j) || fCalibration->GetRingMap(j) > highRing) continue; // Test if the single energy under consideration is with no bigger // then 3 variances from the mean. Note that we do a negated // test. if (rdo->GetSingleEnergy(j) > energyAvg + 3 * energyVar) continue; // increase the number of hits, the energy and the // multiplicity. if(rdo->GetSingleMultiplicity(j)>0) { accHits++; accEnergyAvg += rdo->GetSingleEnergy(j); accMultAvg += rdo->GetSingleMultiplicity(j); } } // Find the average accEnergyAvg = 0; accMultAvg = 0; if(accHits>0) { accEnergyAvg /= accHits; accMultAvg /= accHits; } // cout << "Acc <E>: " << accEnergyAvg // << " Acc <N>: " << accMultAvg << endl; // Now calculate the ring multiplicity and energy Int_t singleRingHits = 0; Double_t ringFinalEnergy = 0; Double_t ringFinalMult = 0; for (k = 0; k < nSingles; k++) { if (fCurrentRing != fCalibration->GetRingMap(k)) continue; // A single is ignored (rather subsituted with a value) if it's // more then 3 variances away from the mean over the ring and // it's immediate neighbors. if (rdo->GetSingleEnergy(k) > energyAvg + 3 * energyVar) continue; if(rdo->GetSingleMultiplicity(k)>0) { singleRingHits++; ringFinalEnergy += rdo->GetSingleEnergy(k); ringFinalMult += rdo->GetSingleMultiplicity(k); } } // Pad with the mean over the ring for those singles not accepted. rdo->SetRingEnergy(fCurrentRing, TMath::Max(0., ringFinalEnergy + accEnergyAvg * (nSingles / nRings - singleRingHits))); rdo->SetRingMultiplicity(fCurrentRing, TMath::Max(0., ringFinalMult + accMultAvg * (nSingles / nRings - singleRingHits))); rdo->SetArrayEnergy(rdo->GetArrayEnergy() + rdo->GetRingEnergy(fCurrentRing)); rdo->SetArrayMultiplicity(rdo->GetArrayMultiplicity() + rdo->GetRingMultiplicity(fCurrentRing)); } } //__________________________________________________________________ void BrMultRdoModule::CalibrateRingsFractionArray(BrMultRdo* rdo) { // PRIVATE METHOD: // Correct for outliers. // The "Fraction Array" method tags single that have a calibrated ADC // signal greater then 10% of the total array calibrated ADC signal // as an outlier. The outliers energy signal is substituted with the // mean over the rest of the array. // Loop variables, must be decalered here so that MSVC++/KCC // doesn't get confused - poor sods, they should stick to ANSI/ISO // standard, alias they don't! Int_t i = 0; Int_t j = 0; // Get some parameters on this detector Int_t nRings = fParameters->GetNoRings(); Int_t nSingles = fParameters->GetNoModules(); for (i = 0; i < nRings; i++) { fCurrentRing = i; // Now calculate the ring multiplicity and energy Int_t singleRingHits = 0; Double_t ringFinalEnergy = 0; Double_t ringFinalMult = 0; for (j = 0; j < nSingles; j++) { if (fCurrentRing != fCalibration->GetRingMap(j)) continue; // A single is ignored it contributes more then 10% to the sum of // the calibrated ADC for the full array. if (rdo->GetSingleCalibratedAdc(j) > .1 * rdo->GetArrayCalibratedAdc()) continue; singleRingHits++; ringFinalEnergy += rdo->GetSingleEnergy(j); if(rdo->GetSingleMultiplicity(j)>0) ringFinalMult += rdo->GetSingleMultiplicity(j); } // Multiply with <possible hits>/<actual hits>, to correct for // outliers. rdo->SetRingEnergy(fCurrentRing, TMath::Max(0., ringFinalEnergy * nSingles / nRings / singleRingHits)); rdo->SetRingMultiplicity(fCurrentRing, TMath::Max(0., ringFinalMult * nSingles / nRings / singleRingHits)); rdo->SetArrayEnergy(rdo->GetArrayEnergy() + rdo->GetRingEnergy(fCurrentRing)); rdo->SetArrayMultiplicity(rdo->GetArrayMultiplicity() + rdo->GetRingMultiplicity(fCurrentRing)); } } //__________________________________________________________________ void BrMultRdoModule::CalibrateArray(BrMultRdo* rdo) { // PRIVATE METHOD: // Make the calibrated data for the full detector. //TESTING if (DebugLevel() > 25) cout << GetName() << " Before correction we have array mult " << rdo->GetArrayMultiplicity() << flush; rdo->SetArrayMultiplicity(TMath::Max(0., rdo->GetArrayMultiplicity() * CalibrateArrayMultiplicity())); if (DebugLevel() > 25) cout << " after " << rdo->GetArrayMultiplicity() << endl; // Fill histograms if (HistOn()) { fArrayEnergyHisto->Fill(rdo->GetArrayEnergy()); fArrayMultHisto->Fill(rdo->GetArrayMultiplicity()); fVtxVsMultHisto->Fill(fCurrentVtxZ, rdo->GetArrayMultiplicity()); } } //__________________________________________________________________ Double_t BrMultRdoModule::CalibrateEta() { // PRIVATE METHOD: // Calculate the pseudo-rapidity for the current ring. if (fCurrentVtxZ == FLT_MAX) { if (DebugLevel() > 25) Warning("CalibrateEta()", "bad vertex"); return FLT_MAX; } Double_t distZ = fCurrentRingZ - fCurrentVtxZ; if (TMath::Abs(distZ) < 1e-6) // distance between center of single and vertex is almost zero, so // we assume it is zero, and set the psuedo rapidity to ZERO return 0; Double_t theta = TMath::ATan2(fCurrentRingR, distZ); return - TMath::Log(TMath::Tan(theta/2)); } //____________________________________________________________________ Double_t BrMultRdoModule::CalibrateAvgEta(Double_t & ringScale) { //PRIVATE METHOD // Calculate average pseudorapidity accounting for depth of detectors. // The method returns the geometric scaling factor (ringScale), where // dN/deta = ringScale * (# of primary hits). // // |--- width -------| // // ----------------- --- // | | | // | si or tile | Depth // | | | // ----------------- --- --- // | // -- delta -eta | // | // -- front // | // -- | // | // ------------------XX (vertex) -------------------- // // Assuming the azmulthal acceptance of a given detector is delta-phi, // the geometric factor to obtain dN/dEta for that detector is: // // ringScale = ( 2*pi/delta-phi ) * 1/delta-eta // if (fCurrentVtxZ == FLT_MAX) { if (DebugLevel() > 25) Warning("CalibrateAvgEta()", "bad vertex"); return FLT_MAX; } ringScale = 0; Double_t x,z; Double_t thmin, thmax, thavg; Double_t etamin,etamax,etaav; Double_t front; Double_t back; Double_t width = fCalibration->GetRingDetectorWidth(); Double_t height = fCalibration->GetRingDetectorHeight(); z = fCurrentRingZ - width/2 - fCurrentVtxZ; front = fCalibration->GetRingDetectorFront(); back = front + fCalibration->GetRingDetectorDepth(); x = front; if(z > 0 ) x = back; thmin = TMath::ATan2(x,z); z = fCurrentRingZ + width/2 - fCurrentVtxZ; x = front; if(z < 0 ) x = back; thmax = TMath::ATan2(x,z); etamin = TMath::Tan(thmin/2.0); if(etamin <= 0) { etamin = 10.0; return etamin; } else etamin = - TMath::Log(etamin); etamax = TMath::Tan(thmax/2.0); if(etamax <= 0) { etamax = 10.0; return etamax; } else etamax = - TMath::Log(etamax); ringScale = 3.1415926/TMath::ATan2( height/2., front); ringScale = ringScale/TMath::Abs(etamin-etamax); if(ringScale<0||ringScale>1000) ringScale = 0; etaav = (etamin+etamax)/2; if(fabs(etaav)>6) etaav = -10; return etaav; } //__________________________________________________________________ Double_t BrMultRdoModule::CalibrateAdc(Int_t adc, Bool_t& below) { // PRIVATE METHOD: // Calibrate the ADC from DAQ. That is, do the pedestal subtraction // if above threshold. below = kFALSE; if (fCurrentRing < 0) { if (DebugLevel() > 25) Warning("CalibrateAdc","A null ring for single"); return -1; } if (fCurrentSingle < 0) { if (DebugLevel() > 25) Warning("CalibrateAdc","A null single %d", fCurrentSingle); return -1; } Double_t pedestal = fCalibration->GetPedestal(fCurrentSingle); Int_t gap = fCalibration->GetAdcGap(fCurrentSingle) ; Double_t threshold = fThresholdFactor * fCalibration->GetPedestalWidth(fCurrentSingle) + pedestal; // Test if this single is active, and signal is above // cout << "cholm thres: " << threshold << endl; if (adc < threshold) { if (Verbose() > 25) cout << "Adc " << fCurrentSingle << " below threshold: " << adc << " < " << threshold << " (=" << fThresholdFactor << " * " << fCalibration->GetPedestalWidth(fCurrentSingle) << " + " << pedestal << ")" << endl; below = kTRUE; } // Correct for ADC gab if (fAdcGapLimit > 0 && adc > fAdcGapLimit) adc -= gap + 1; // subtract pedestal return Double_t(adc - pedestal); } //__________________________________________________________________ Double_t BrMultRdoModule::CalibrateEnergy(Double_t calAdc) { // PRIVATE METHOD: // The energy in a given single is calculated as // // sin(theta + tilt) // E = calADC * gain * ----------------- (1) // sin(theta) // // The third factor in (1) is a correction for the tilt of the // detector. // Given that // // sin(theta + tilt) = sin(theta)*cos(tilt) + sin(tilt)cos(theta) // // and that // // ringR // sin(theta) = ----------------------- // sqrt(distZ^2 + ringR^2) // // distZ // cos(theta) = ----------------------- // sqrt(distZ^2 + ringR^2) // // distZ // cotan(theta) = ----- // ringR // // // ------- ---+--- ------- ------- // /| // / | ringR // / | // =======*===+================== // distZ // // where distZ = ringZ - vtxZ, we get // // sin(theta - tilt) sin(theta)*cos(tilt) - sin(tilt)cos(theta) // ----------------- = ------------------------------------------ // sin(theta) sin(theta) // // cos(theta) // = cos(tilt) - ---------- * sin(tilt) // sin(theta) // // = cos(tilt) - cotan(theta) * sin(tilt) // // distZ // = cos(tilt) - ----- * sin(tilt) // ringR // // so we get // // // distZ // E = calADC * gain * ( cos(tilt) - ----- * sin(tilt) ) // ringR // // We still need to correct for the vertex position, but that is // done in the Event method. // if (fCurrentRing < 0) { if (DebugLevel() > 25) Warning("CalibrateEnergy","A null ring for single"); return 0; } if (fCurrentSingle < 0) { if (DebugLevel() > 25) Warning("CalibrateEnergy","A null single %d", fCurrentSingle); return 0; } if (fCurrentVtxZ == FLT_MAX) { if (DebugLevel() > 25) Warning("CalibrateEnergy","A null single vertex"); return 0; } Double_t tilt = fCalibration->GetTilt(fCurrentSingle); Double_t gain = fCalibration->GetAdcGain(fCurrentSingle); Double_t distZ = fCurrentRingZ - fCurrentVtxZ; // Calculate the energy return calAdc * gain * (TMath::Cos(tilt) - distZ / fCurrentRingR * TMath::Sin(tilt)) / 1000.; } //__________________________________________________________________ Double_t BrMultRdoModule::CalibrateMultiplicity(Double_t ringEta) { // PRIVATE METHOD: // Calculate the Multiplicity from the energy deposited in a single if (ringEta == FLT_MAX) { if (DebugLevel() > 10) Warning("CalibrateMultiplicity", "bad eta value %f", ringEta); return 0; } Int_t order = fCalibration->GetConversionFuncOrder(); if (order < 0) { if (DebugLevel() > 10) Warning("CalibrateMultiplicity","A null function"); return 0; } // Calculate the conversion factor Double_t mip = fCalibration->GetConversionFuncPar(fCurrentRing, 0); for (Int_t i = 1; i < order + 1; i++) mip += fCalibration->GetConversionFuncPar(fCurrentRing, i) * TMath::Power(ringEta, i); return mip; } //____________________________________________________________________ // // $Log: BrMultRdoModule.cxx,v $ // Revision 1.20 2002/05/08 16:30:55 cholm // Added 3 histograms, one of the sum of the multiplicity in each eta bin, // one with the number of fills, and a third, the ratio of these two. Unless // something tricky is going on, this should be our dN/deta distribution for // the TMA and SMA. Also fixed a few warning (probably GCC3 errors) in the // Si and Tile modules. // // Revision 1.19 2002/04/19 21:33:15 hito // Minor change for compiling with gcc3. // // Revision 1.18 2002/03/05 20:28:49 sanders // corrected dndeta calculation // // Revision 1.17 2002/02/19 21:26:39 sanders // modified pointer argument to call by reference // // Revision 1.16 2002/01/27 16:12:52 cholm // Use pass-by-reference ('&') rather than pointers ('*') in the method // BrMultRdoModule::CalibrateAdc, since that's the proper C++ way to go // about these things. // // Revision 1.15 2002/01/26 16:17:33 sanders // Modified point where the energy threshold is applied. Delaying the // threshold cut until after the energy calibration is necessary for // checking the energy calibration. // Also corrected an indexing error that was cutting off the first adc // channel in the histograms. // // Revision 1.14 2002/01/04 21:30:11 sanders // Added "Setters" for range of raw adc values. // // Revision 1.13 2001/12/07 18:12:08 cholm // Better range on ADC Y axis // // Revision 1.12 2001/12/07 14:34:54 sanders // another check added for spurious values being added into summed mult. // // Revision 1.11 2001/12/06 23:38:29 sanders // placed check to avoid summing in (-1) that don't exist // when finding sum multiplicity. // // Revision 1.10 2001/12/06 00:35:08 sanders // Revised method for obtaining average pseudorapidity and // geometric scaling factors for rings. // // Revision 1.9 2001/11/25 16:00:28 sanders // fixed missing semicolon... // // Revision 1.8 2001/11/25 15:52:05 sanders // Removed gap closing factor that was introduced by mistake. // Added checks for vertex being in range for dNdEta and AvgEta calculations. // // Revision 1.7 2001/11/24 15:47:08 sanders // Added code to calculate dNdEta values for rings (average over // individual dNdEta values calculated for each detector) and // HIJING weighted pseudorapidity. // // Revision 1.6 2001/11/19 04:07:48 sanders // Revised rdo modules for 2001 si and tile calibrations. // // Revision 1.5 2001/11/02 15:14:51 cholm // Changed the eta range to [-3,3] // // Revision 1.4 2001/10/29 21:01:08 cholm // Corrected a very bad mistake in the previous commit: A manager was // initialised in the module - BAD. Managers may not be initialised in // modules, only thier services may be used. // // Revision 1.3 2001/09/20 13:44:10 cholm // Removed the use of temp ASCII files, since not needed anymore // // Revision 1.2 2001/08/03 21:36:30 zdc // The methods to read low-gain ZDC ADC have been added to BrZdcRdoModule // // Revision 1.1.1.1 2001/06/21 14:55:08 hagel // Initial revision of brat2 // // Revision 1.3 2001/06/06 19:39:31 sanders // Tilt angle correction for individual detector elements fixed. The // <si/tile> centrality function parameters changed accordingly. // // Revision 1.2 2001/06/01 15:49:18 cholm // Fix of some default values from -1 to zero (gave bad results) // BrSiCalibration had a bug in returning the correction function. // Tile|Rdo modules had some minor bugs, and I did some improvments // Mult Rdo module had some serious bugs, but should be ok now. // BrMultRdo is fixed a lot too. // // Revision 1.1 2001/05/31 01:46:08 cholm // Second rewamp of this directory. All RDO modules use the common // BrMultRdoModule, and they both write BrMultRdo Objects. Also introduced // ABC for Br[Tile|Si][Parameters|Calibration] since they have a lot in common, // and the code can be made more efficient this way. // // A possible further thing to do, is to make an ABC for the CentModules, and // the corresponding calibration modules, since they are very VERY similar. // However, the current module BrMultCentModule is in the way. Need to talk // to Steve about that. // // |
||||||
|
This page automatically generated by script docBrat by Christian Holm |
Copyright ; 2002 BRAHMS Collaboration
<brahmlib@rcf.rhic.bnl.gov>
|