// fSignals->GetXXX..()
// etc.
//________________________________________________________________________
+#include <string>
+#include <sstream>
+#include <fstream>
#include "TProfile.h"
#include "TFile.h"
#include "AliRawReader.h"
-#include "AliRawEventHeaderBase.h"
#include "AliCaloRawStreamV3.h"
//The include file
fMapping(NULL),
fRunNumber(-1),
fStartTime(0),
- fAmpCut(50),
- fReqFractionAboveAmpCutVal(0.8),
+ fAmpCut(40), // min. 40 ADC counts as default
+ fReqFractionAboveAmpCutVal(0.6), // 60% in a strip, per default
fReqFractionAboveAmp(kTRUE),
fHour(0),
fLatestHour(0),
//_____________________________________________________________________
void AliCaloCalibSignal::ResetInfo()
-{
+{ // reset trees and counters
Zero(); // set all counters to 0
DeleteTrees(); // delete previous stuff
CreateTrees(); // and create some new ones
fNEvents = 0;
fNAcceptedEvents = 0;
- // Set the number of points for each tower: Amp vs. Time
+ // Set the number of points for each tower: Amp vs. Time
memset(fNHighGain, 0, sizeof(fNHighGain));
memset(fNLowGain, 0, sizeof(fNLowGain));
// and LED reference
}
//_____________________________________________________________________
-Bool_t AliCaloCalibSignal::CheckFractionAboveAmp(int *AmpVal, int nTotChan)
-{
- int nAbove = 0;
+Bool_t AliCaloCalibSignal::CheckFractionAboveAmp(const int *iAmpVal,
+ int resultArray[])
+{ // check fraction of towers, per column, that are above amplitude cut
+ Bool_t returnCode = false;
- int TowerNum = 0;
+ int iTowerNum = 0;
+ double fraction = 0;
for (int i = 0; i<fModules; i++) {
for (int j = 0; j<fColumns; j++) {
+ int nAbove = 0;
for (int k = 0; k<fRows; k++) {
- TowerNum = GetTowerNum(i,j,k);
- if (AmpVal[TowerNum] > fAmpCut) {
+ iTowerNum = GetTowerNum(i,j,k);
+ if (iAmpVal[iTowerNum] > fAmpCut) {
nAbove++;
}
}
+ resultArray[i*fColumns +j] = 0; // init. to denied
+ if (nAbove > 0) {
+ fraction = (1.0*nAbove) / fRows;
+ /*
+ printf("DS mod %d col %d nAbove %d fraction %3.2f\n",
+ i, j, nAbove, fraction);
+ */
+ if (fraction > fReqFractionAboveAmpCutVal) {
+ resultArray[i*fColumns + j] = nAbove;
+ returnCode = true;
+ }
+ }
}
- }
-
- double fraction = (1.0*nAbove) / nTotChan;
+ } // modules loop
- if (fraction > fReqFractionAboveAmpCutVal) {
- return true;
+ return returnCode;
+}
+
+// Parameter/cut handling
+//_____________________________________________________________________
+void AliCaloCalibSignal::SetParametersFromFile(const char *parameterFile)
+{ // set parameters from file
+ static const string delimitor("::");
+
+ // open, check input file
+ ifstream in( parameterFile );
+ if( !in ) {
+ printf("in AliCaloCalibSignal::SetParametersFromFile - Using default/run_time parameters.\n");
+ return;
+ }
+
+ // Note: this method is a bit more complicated than it really has to be
+ // - allowing for multiple entries per line, arbitrary order of the
+ // different variables etc. But I wanted to try and do this in as
+ // correct a C++ way as I could (as an exercise).
+
+ // read in
+ char readline[1024];
+ while ((in.rdstate() & ios::failbit) == 0 ) {
+
+ // Read into the raw char array and then construct a string
+ // to do the searching
+ in.getline(readline, 1024);
+ istringstream s(readline);
+
+ while ( ( s.rdstate() & ios::failbit ) == 0 ) {
+
+ string keyValue;
+ s >> keyValue;
+
+ // check stream status
+ if( s.rdstate() & ios::failbit ) break;
+
+ // skip rest of line if comments found
+ if( keyValue.substr( 0, 2 ) == "//" ) break;
+
+ // look for "::" in keyValue pair
+ size_t position = keyValue.find( delimitor );
+ if( position == string::npos ) {
+ printf("wrong format for key::value pair: %s\n", keyValue.c_str());
+ }
+
+ // split keyValue pair
+ string key( keyValue.substr( 0, position ) );
+ string value( keyValue.substr( position+delimitor.size(),
+ keyValue.size()-delimitor.size() ) );
+
+ // check value does not contain a new delimitor
+ if( value.find( delimitor ) != string::npos ) {
+ printf("wrong format for key::value pair: %s\n", keyValue.c_str());
+ }
+
+ // debug: check key value pair
+ // printf("AliCaloCalibSignal::SetParametersFromFile - key %s value %s\n", key.c_str(), value.c_str());
+
+ // if the key matches with something we expect, we assign the new value
+ if ( (key == "fAmpCut") || (key == "fReqFractionAboveAmpCutVal") || (key == "fSecInAverage") ) {
+ istringstream iss(value);
+ printf("AliCaloCalibSignal::SetParametersFromFile - key %s value %s\n", key.c_str(), value.c_str());
+
+ if (key == "fAmpCut") {
+ iss >> fAmpCut;
+ }
+ else if (key == "fReqFractionAboveAmpCutVal") {
+ iss >> fReqFractionAboveAmpCutVal;
+ }
+ else if (key == "fSecInAverage") {
+ iss >> fSecInAverage;
+ }
+ } // some match found/expected
+
+ }
}
- else return false;
+
+ in.close();
+ return;
+
+}
+
+//_____________________________________________________________________
+void AliCaloCalibSignal::WriteParametersToFile(const char *parameterFile)
+{ // write parameters to file
+ static const string delimitor("::");
+ ofstream out( parameterFile );
+ out << "// " << parameterFile << endl;
+ out << "fAmpCut" << "::" << fAmpCut << endl;
+ out << "fReqFractionAboveAmpCutVal" << "::" << fReqFractionAboveAmpCutVal << endl;
+ out << "fSecInAverage" << "::" << fSecInAverage << endl;
+
+ out.close();
+ return;
}
//_____________________________________________________________________
// if fMapping is NULL the rawstream will crate its own mapping
AliCaloRawStreamV3 rawStream(rawReader, fCaloString, (AliAltroMapping**)fMapping);
- return ProcessEvent( &rawStream, (AliRawEventHeaderBase*)rawReader->GetEventHeader() );
+ return ProcessEvent( &rawStream, rawReader->GetTimestamp() );
}
//_____________________________________________________________________
-Bool_t AliCaloCalibSignal::ProcessEvent(AliCaloRawStreamV3 *in, AliRawEventHeaderBase *aliHeader)
+Bool_t AliCaloCalibSignal::ProcessEvent(AliCaloRawStreamV3 *in, UInt_t Timestamp)
{
// Method to process=analyze one event in the data stream
if (!in) return kFALSE; //Return right away if there's a null pointer
fNEvents++; // one more event
// use maximum numbers to set array sizes
- int AmpValHighGain[fgkMaxTowers];
- int AmpValLowGain[fgkMaxTowers];
- memset(AmpValHighGain, 0, sizeof(AmpValHighGain));
- memset(AmpValLowGain, 0, sizeof(AmpValLowGain));
+ int iAmpValHighGain[fgkMaxTowers];
+ int iAmpValLowGain[fgkMaxTowers];
+ memset(iAmpValHighGain, 0, sizeof(iAmpValHighGain));
+ memset(iAmpValLowGain, 0, sizeof(iAmpValLowGain));
// also for LED reference
- int LEDAmpVal[fgkMaxRefs * 2]; // factor 2 is for the two gain values
- memset(LEDAmpVal, 0, sizeof(LEDAmpVal));
+ int iLEDAmpVal[fgkMaxRefs * 2]; // factor 2 is for the two gain values
+ memset(iLEDAmpVal, 0, sizeof(iLEDAmpVal));
int sample; // temporary value
int gain = 0; // high or low gain
int nLowChan = 0;
int nHighChan = 0;
- int TowerNum = 0; // array index for regular towers
- int RefNum = 0; // array index for LED references
+ int iTowerNum = 0; // array index for regular towers
+ int iRefNum = 0; // array index for LED references
// loop first to get the fraction of channels with amplitudes above cut
// counters
int max = AliEMCALGeoParams::fgkSampleMin, min = AliEMCALGeoParams::fgkSampleMax; // min and max sample values
-
+ int nsamples = 0;
+
while (in->NextBunch()) {
const UShort_t *sig = in->GetSignals();
+ nsamples += in->GetBunchLength();
for (Int_t i = 0; i < in->GetBunchLength(); i++) {
sample = sig[i];
} // loop over samples in bunch
} // loop over bunches
+ if (nsamples > 0) { // this check is needed for when we have zero-supp. on, but not sparse readout
+
gain = -1; // init to not valid value
//If we're here then we're done with this tower
if ( in->IsLowGain() ) {
if ( in->IsHighGain() || in->IsLowGain() ) { // regular tower
// get tower number for AmpVal array
- TowerNum = GetTowerNum(arrayPos, in->GetColumn(), in->GetRow());
+ iTowerNum = GetTowerNum(arrayPos, in->GetColumn(), in->GetRow());
if (gain == 0) {
// fill amplitude into the array
- AmpValLowGain[TowerNum] = max - min;
+ iAmpValLowGain[iTowerNum] = max - min;
nLowChan++;
}
else if (gain==1) {//fill the high gain ones
// fill amplitude into the array
- AmpValHighGain[TowerNum] = max - min;
+ iAmpValHighGain[iTowerNum] = max - min;
nHighChan++;
}//end if gain
} // regular tower
else if ( in->IsLEDMonData() ) { // LED ref.;
// strip # is coded is 'column' in the channel maps
- RefNum = GetRefNum(arrayPos, in->GetColumn(), gain);
- LEDAmpVal[RefNum] = max - min;
+ iRefNum = GetRefNum(arrayPos, in->GetColumn(), gain);
+ iLEDAmpVal[iRefNum] = max - min;
} // end of LED ref
-
+
+ } // nsamples>0 check, some data found for this channel; not only trailer/header
} // end while over channel
}//end while over DDL's, of input stream
in->Reset(); // just in case the next customer forgets to check if the stream was reset..
// now check if it was a led event, only use high gain (that should be sufficient)
+
+ // by default all columns are accepted (init check to > 0)
+ int checkResultArray[AliEMCALGeoParams::fgkEMCALModules * AliEMCALGeoParams::fgkEMCALCols];
+ for (int ia=0; ia<(AliEMCALGeoParams::fgkEMCALModules * AliEMCALGeoParams::fgkEMCALCols); ia++) {
+ checkResultArray[ia] = 1;
+ }
+
if (fReqFractionAboveAmp) {
bool ok = false;
if (nHighChan > 0) {
- ok = CheckFractionAboveAmp(AmpValHighGain, nHighChan);
+ ok = CheckFractionAboveAmp(iAmpValHighGain, checkResultArray);
}
if (!ok) return false; // skip event
}
fNAcceptedEvents++; // one more event accepted
if (fStartTime == 0) { // if start-timestamp wasn't set,we'll pick it up from the first event we encounter
- fStartTime = aliHeader->Get("Timestamp");
+ fStartTime = Timestamp;
}
- fHour = (aliHeader->Get("Timestamp")-fStartTime)/(double)fgkNumSecInHr;
+ fHour = (Timestamp - fStartTime)/(double)fgkNumSecInHr;
if (fLatestHour < fHour) {
fLatestHour = fHour;
}
// it is a led event, now fill TTree
for(int i=0; i<fModules; i++){
for(int j=0; j<fColumns; j++){
+ if (checkResultArray[i*fColumns + j] > 0) { // column passed check
for(int k=0; k<fRows; k++){
- TowerNum = GetTowerNum(i, j, k);
+ iTowerNum = GetTowerNum(i, j, k);
- if(AmpValHighGain[TowerNum]) {
- fAmp = AmpValHighGain[TowerNum];
+ if(iAmpValHighGain[iTowerNum]) {
+ fAmp = iAmpValHighGain[iTowerNum];
fChannelNum = GetChannelNum(i,j,k,1);
- fTreeAmpVsTime->Fill();//fChannelNum,fHour,AmpValHighGain[TowerNum]);
- fNHighGain[TowerNum]++;
+ fTreeAmpVsTime->Fill();//fChannelNum,fHour,AmpValHighGain[iTowerNum]);
+ fNHighGain[iTowerNum]++;
}
- if(AmpValLowGain[TowerNum]) {
- fAmp = AmpValLowGain[TowerNum];
+ if(iAmpValLowGain[iTowerNum]) {
+ fAmp = iAmpValLowGain[iTowerNum];
fChannelNum = GetChannelNum(i,j,k,0);
- fTreeAmpVsTime->Fill();//fChannelNum,fHour,AmpValLowGain[TowerNum]);
- fNLowGain[TowerNum]++;
+ fTreeAmpVsTime->Fill();//fChannelNum,fHour,AmpValLowGain[iTowerNum]);
+ fNLowGain[iTowerNum]++;
}
} // rows
+ } // column passed check
} // columns
+ // We also do the activity check for LEDRefs/Strips, but need to translate between column
+ // and strip indices for that; based on these relations:
+ // iStrip = AliEMCALGeoParams::GetStripModule(iSM, iCol);
+ // iStrip = (iSM%2==0) ? iCol/2 : AliEMCALGeoParams::fgkEMCALLEDRefs - 1 - iCol/2;
+ // which leads to
+ // iColFirst = (iSM%2==0) ? iStrip*2 : (AliEMCALGeoParams::fgkEMCALLEDRefs - 1 - iStrip)*2;
+
// also LED refs
for(int j=0; j<fLEDRefs; j++){
+ int iColFirst = (i%2==0) ? j*2 : (AliEMCALGeoParams::fgkEMCALLEDRefs - 1 - j)*2; //CHECKME!!!
+ if ( (checkResultArray[i*fColumns + iColFirst]>0) || (checkResultArray[i*fColumns + iColFirst + 1]>0) ) { // at least one column in strip passed check
for (gain=0; gain<2; gain++) {
fRefNum = GetRefNum(i, j, gain);
- if (LEDAmpVal[fRefNum]) {
- fAmp = LEDAmpVal[fRefNum];
+ if (iLEDAmpVal[fRefNum]) {
+ fAmp = iLEDAmpVal[fRefNum];
fTreeLEDAmpVsTime->Fill();//fRefNum,fHour,fAmp);
fNRef[fRefNum]++;
}
- }
+ } // gain
+ } // at least one column in strip passed check
}
} // modules
memset(profile, 0, sizeof(profile));
char name[200]; // for profile id and title
- int TowerNum = 0;
+ int iTowerNum = 0;
for (int i = 0; i<fModules; i++) {
for (int ic=0; ic<fColumns; ic++){
for (int ir=0; ir<fRows; ir++) {
- TowerNum = GetTowerNum(i, ic, ir);
+ iTowerNum = GetTowerNum(i, ic, ir);
// high gain
- if (fNHighGain[TowerNum] > 0) {
+ if (fNHighGain[iTowerNum] > 0) {
fChannelNum = GetChannelNum(i, ic, ir, 1);
sprintf(name, "profileChan%d", fChannelNum);
profile[fChannelNum] = new TProfile(name, name, numProfBins, timeMin, timeMax, "s");
}
// same for low gain
- if (fNLowGain[TowerNum] > 0) {
+ if (fNLowGain[iTowerNum] > 0) {
fChannelNum = GetChannelNum(i, ic, ir, 0);
sprintf(name, "profileChan%d", fChannelNum);
profile[fChannelNum] = new TProfile(name, name, numProfBins, timeMin, timeMax, "s");
return kTRUE;
}
+
+//_____________________________________________________________________
+Bool_t AliCaloCalibSignal::DecodeChannelNum(const int chanId,
+ int *imod, int *icol, int *irow, int *igain) const
+{ // return the module, column, row, and gain for a given channel number
+ *igain = chanId/(fModules*fColumns*fRows);
+ *imod = (chanId/(fColumns*fRows)) % fModules;
+ *icol = (chanId/fRows) % fColumns;
+ *irow = chanId % fRows;
+ return kTRUE;
+}
+
+//_____________________________________________________________________
+Bool_t AliCaloCalibSignal::DecodeRefNum(const int refId,
+ int *imod, int *istripMod, int *igain) const
+{ // return the module, stripModule, and gain for a given reference number
+ *igain = refId/(fModules*fLEDRefs);
+ *imod = (refId/(fLEDRefs)) % fModules;
+ *istripMod = refId % fLEDRefs;
+ return kTRUE;
+}
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff