//*-- Authors: Leszek Kidon & Jacek Otwinowski //*-- Last Modified: 08/07/2001 (Marcin Jaskula) //*-- Last Modified: 17/11/2005 (Jacek Otwinowski) //*-- Last Modified: 05/12/2006 (Jacek Otwinowski) //_HADES_CLASS_DESCRIPTION ///////////////////////////////////////////////////////////////////////////// // // HShowerPadDigitizer // // This class digitizes the shower pads. For each fired wire it calculates all // the pads which the charge couples to. For each pad the track numbers of the // particles that fired the pad are stored in a linear category (HShowerTrack). // This category is sortable, in particular the track numbers are sorted by // the respective pad address. // All the fired pads are stored in a matrix category (catShowerRawMatr, // this matrix category is used only for the simulated data). // // The Shower digitization is split into several tasks as shown // in the flow diagram below. // // ---------------------- // // | HShowerUnpacker | // // | (embedding mode) | \ // // | | \ ------------------ // // ---------------------- | | HGeantShower | // // | ------------------\ // // | \ // // | ------------------ \-------------> ---------------------- // // | | HGeantWire | <------------ | HShowerHitDigitizer | // // | ------------------\ ---------------------- // // | \ // // ------------- ------------------ \-------------> ----------------------- // // -- | HShowerRaw | | HShowerRawMatr | <------------ | HShowerPadDigitizer | // // | ------------- ------------------\ |( creates track objects| // // | \ | for real tracks in | // // ---------------------- ------------------ \ | embedding mode too) | // // | HShowerCalibrater | | HShowerTrack | \<------------ ----------------------- // // | (embedding mode) | ------------------\ \ // // ---------------------- \ \ ----------------------- // // | ------------------ \ ----------> | HShowerCopy | // // -------------------> | HShowerCal | \<------------ |(add charge of real hit| // // ------------------\ \ | in embedding too ) | // // \ \ ----------------------- // // ------------------ ----\---------> ----------------------- // // | HShowerHitHdr | <--\--------- | HShowerHitFinder | // // ------------------ \ ----------------------- // // ------------------ \ | // // | HShowerPID | <-----\-------| // // ------------------ \ | // // ------------------ \ | // // | HShowerHit | <--------\----| // // ------------------ < \ // // \ \ // // \-------------> ------------------------ // // | HShowerHitTrackMatcher | // // ------------------------ // // // // // // In the case of TRACK EMBEDDING of simulated tracks into // experimental data the real data are written by the HShowerUnpacker into // HShowerRaw category. The real hits are taken into // account by the digitizer (adding of charges). The embedding mode is recognized // automatically by analyzing the // gHades->getEmbeddingMode() flag. // Mode ==0 means no embedding // ==1 realistic embedding (first real or sim hit makes the game) // ==2 keep GEANT tracks (as 1, but GEANT track numbers will always // win against real data. besides the tracknumber the output will // be the same as in 1) // //------------------------------------------------- // Used to suppress the secondaries created in the // SHOWER itself. // 0 = realistic (secondaries included) // 1 primary particle is stored // 2 the first track number entering the SHOWER in SAME SECTOR is stored // 3 the first track number entering the TOFINO in SAME SECTOR is stored // or the primary track if no TOFINO was found // 4 (default) the first track number entering the TOFINO in SAME SECTOR is stored // or the first track in SHOWER if no TOFINO was found // // // The mode can be selected by static void HShowerPadDigitizer::setModeTrack(Int_t mode) // ////////////////////////////////////////////////////////////////////////////// #include "hshowerpaddigitizer.h" #include #include #include #include "hevent.h" #include "hspectrometer.h" #include "hdetector.h" #include "hshowerdetector.h" #include "hcategory.h" #include "hmatrixcategory.h" #include "hlinearcategory.h" #include "hlocation.h" #include "hshowergeantwire.h" #include "hshowercal.h" #include "hshowertrack.h" #include "hshowerdigipar.h" #include "hshowergeometry.h" #include "hdebug.h" #include "hades.h" #include "hgeantshower.h" #include "hgeanttof.h" #include "hgeantrpc.h" #include "hgeantkine.h" #include "showerdef.h" #include "TMath.h" #include #define MIN(A,B) (((A) <= (B)) ? (A) : (B)) #define MAX(A,B) (((A) >= (B)) ? (A) : (B)) #define TABLE_SIZE 10 ClassImp(HShowerPadDigitizer) Int_t HShowerPadDigitizer::modeTrack = 4; HShowerPadDigitizer::HShowerPadDigitizer(const Text_t *name,const Text_t *title) : HShowerDigitizer(name,title) { fTrackIter = NULL; fChannelCoeff = 256.0 / 60.0; // QDC calibration factor fShowerRawMatrIter = NULL; isRpc = kFALSE; } HShowerPadDigitizer::HShowerPadDigitizer() { fTrackIter = NULL; fShowerRawMatrIter = NULL; isRpc = kFALSE; } HShowerPadDigitizer::~HShowerPadDigitizer(void) { if (fTrackIter) delete fTrackIter; if (fShowerRawMatrIter) delete fShowerRawMatrIter; } Bool_t HShowerPadDigitizer::init() { // creates ShowerGeantWire(input), ShowerRawMatr(output) and ShowerTrack // categories and adds them to current event // creates an iterator for ShowerTrack category printf("initialization of shower pad digitizer \n"); HCategory *pCat; HShowerDetector *pShowerDet = (HShowerDetector*)gHades->getSetup()->getDetector("Shower"); pCat=gHades->getCurrentEvent()->getCategory(catShowerGeantWire); if (!pCat) { pCat=pShowerDet->buildCategory(catShowerGeantWire); if (!pCat) return kFALSE; else gHades->getCurrentEvent()->addCategory(catShowerGeantWire, pCat, "Shower"); } setInCat(pCat); pCat=gHades->getCurrentEvent()->getCategory(catShowerRawMatr); if (!pCat) { pCat=pShowerDet->buildCategory(catShowerRawMatr); if (!pCat) return kFALSE; else { gHades->getCurrentEvent()->addCategory(catShowerRawMatr, pCat, "Shower"); } } setOutCat(pCat); if(pCat) fShowerRawMatrIter=(HIterator*)pCat->MakeIterator("native"); if(gHades->getEmbeddingMode()>0) { // we need to get the HShowerCal catgeory // in addition pCat=gHades->getCurrentEvent()->getCategory(catShowerCal); if (!pCat) { Error("init()","No catShowerCal in input ... needed for embedding!"); return kFALSE; } if(pCat) fShowerCalIter= (HIterator*)pCat->MakeIterator("native"); } pCat=gHades->getCurrentEvent()->getCategory(catShowerTrack); if (!pCat) { pCat=pShowerDet->buildCategory(catShowerTrack); if (!pCat) return kFALSE; else gHades->getCurrentEvent()->addCategory(catShowerTrack, pCat, "Shower"); } setTrackCat(pCat); fTrackIter = (HIterator*)getTrackCat()->MakeIterator("native"); //------------------------------------------ // getting sim categories fGeantKineCat = (HLinearCategory*)gHades->getCurrentEvent()->getCategory(catGeantKine); if(!fGeantKineCat){ Error("init()","Received Zero pointer for catGeantKine!"); return kFALSE; } fGeantShowerCat = (HLinearCategory*)gHades->getCurrentEvent()->getCategory(catShowerGeantRaw); if(!fGeantShowerCat){ Error("init()","Received Zero pointer for catShowerGeantRaw!"); return kFALSE; } fGeantTofCat = (HLinearCategory*)gHades->getCurrentEvent()->getCategory(catTofGeantRaw); if(!fGeantTofCat){ Warning("init()","Received Zero pointer for catTofGeantRaw!"); } fGeantRpcCat = (HLinearCategory*)gHades->getCurrentEvent()->getCategory(catRpcGeantRaw); if(!fGeantRpcCat){ Warning("init()","Received Zero pointer for catRpcGeantRaw!"); } if (!fGeantRpcCat && !fGeantTofCat) return kFALSE; if (fGeantRpcCat) isRpc = kTRUE; //------------------------------------------ return HShowerDigitizer::init(); } HShowerPadDigitizer& HShowerPadDigitizer::operator=(HShowerPadDigitizer &c) { // It should have been done return c; } Int_t HShowerPadDigitizer::execute() { trackMap.clear(); // filled in digitize()->digiPads()->updatePad() if(gHades->getEmbeddingMode() > 0) { //--------------------------------------------------------------- // loop over ShowerCal category and create the // corresponding HShowerTrack for the real data objects. HLocation loc; loc.set(4,0,0,0); HShowerCal *pCal; fShowerCalIter->Reset(); HShowerRaw pRaw; while ((pCal = (HShowerCal *)fShowerCalIter->Next()) != 0) { loc[0] = pCal->getSector(); loc[1] = pCal->getModule(); loc[2] = pCal->getRow(); loc[3] = pCal->getCol(); pRaw.setSector(loc[0]); pRaw.setModule(loc[1]); pRaw.setRow(loc[2]); pRaw.setCol(loc[3]); //------------------------------------------------------------ // add pads to track map Int_t addTr = pCal->getAddress(); vector v; v.push_back(gHades->getEmbeddingRealTrackId()); trackMap.insert( make_pair( addTr, v ) ); //------------------------------------------------------------ } } //--------------------------------------------------------------- //--------------------------------------------------------------- TObject *pHit; lNrEvent++; fIter->Reset(); while((pHit = fIter->Next())) { digitize(pHit); } //--------------------------------------------------------------- //--------------------------------------------------------------- // create HShowerTrack object for each unique track of the fired pad HShowerTrack* pShowerTrack; HLocation loc; Int_t s,m,row,col; vector check; for( map >::iterator iter = trackMap.begin(); iter != trackMap.end(); ++iter) { vector& list = iter ->second; Int_t nAddress = iter ->first; s = nAddress /100000; m = (nAddress - 100000*s) / 10000; row = (nAddress - 100000*s - 10000*m) / 100; col = (nAddress - 100000*s - 10000*m - 100*row); loc.set(4,s,m,row,col); Int_t n = list.size(); std::sort( list.begin(), list.end() ); // sort increasing order : real,primary,secondary with TOF,secondary with no TOF check.clear(); for( Int_t i = 0; i < n; ++ i ) { Int_t nTrack = list[i]; if( nTrack > 100000000) { nTrack -= 100000000; } // secondary with no TOF else if(nTrack > 10000) { nTrack -= 10000; } // secondary with TOF else if(nTrack > 0) { } // primary else if(nTrack < 0) { } // real if (i > 0 && find( check.begin(), check.end(), nTrack ) != check.end() ){ // track exists already .... nothing todo } else { // new track ... create HShowerTrack check.push_back(nTrack); pShowerTrack = (HShowerTrack*)((HLinearCategory*)getTrackCat())->getNewSlot(loc); if (pShowerTrack) { pShowerTrack = new(pShowerTrack) HShowerTrack; pShowerTrack->setAddress(nAddress); pShowerTrack->setTrack(nTrack); } } } } //--------------------------------------------------------------- //--------------------------------------------------------------- // gain and efficiency Double_t fQ; HShowerRawMatr *pRawMatr; fShowerRawMatrIter->Reset(); while((pRawMatr = (HShowerRawMatr*)fShowerRawMatrIter->Next())) { // calculate gain fQ = gainCharge(pRawMatr); pRawMatr->setCharge(fQ); // effciency model good for post1 and post2 if(!checkEfficiency(pRawMatr)) { pRawMatr->setCharge(0.0); } } //--------------------------------------------------------------- sort(); // this is used to call the track sort() in pad digitizer return 0; } Bool_t HShowerPadDigitizer::digitize(TObject *pHit) { #if DEBUG_LEVEL>2 gDebuger->enterFunc("HShowerPadDigitizer::execute"); gDebuger->message("Hit cat points to %p",pHit); #endif HShowerGeantWire *pGeantWire = (HShowerGeantWire*) pHit; if (pGeantWire) { digiPads(pGeantWire); } #if DEBUG_LEVEL>2 gDebuger->leaveFunc("HShowerPadDigitizer::execute"); #endif return (pGeantWire) ? kTRUE : kFALSE; } Int_t HShowerPadDigitizer::sort() { // if (((HLinearCategory*)getTrackCat())->IsSortable()) ((HLinearCategory*)getTrackCat())->sort(); return 0; } Int_t HShowerPadDigitizer::findFirstHitInShower(Int_t trackID) { //------------------------------------------------- // Used to suppress the secondaries created in the // SHOWER itself. // 0 = realistic (secondaries included) // 1 primary particle is stored // 2 the first track number entering the SHOWER in SAME SECTOR is stored // 3 the first track number entering the TOFINO in SAME SECTOR is stored // or the primary track if no TOFINO was found // 4 (default) the first track number entering the TOFINO in SAME SECTOR is stored // or the first track in SHOWER if no TOFINO was found Int_t numTrack = trackID; if(numTrack <= 0) return numTrack; // nothing to do for negative track numbers //-------------------------------------- // case 0 : in = out if(modeTrack == 0) { return numTrack; } HGeantShower *poldShower; Int_t first = 0; Int_t parent= 0; HGeantKine* kine = (HGeantKine*)fGeantKineCat->getObject(numTrack - 1); if(kine){ parent = kine->getParentTrack(); if( parent == 0) { return numTrack; } // nothing todo first = kine->getFirstShowerHit(); if(first != -1){ poldShower = (HGeantShower*)fGeantShowerCat->getObject(first); } else { Error("findFirstHitInShower()","No first shower hit!"); return numTrack; } } else { Error("findFirstHitInShower()","Received Zero pointer for kine!"); return numTrack; } if(numTrack != poldShower->getTrack()){ Error("findFirstHitInShower()","First shower hit not same trackID!"); return numTrack; } //-------------------------------------------------------- // return the track number for // the selected option modeTrack Int_t s = poldShower->getSector(); //-------------------------------------- // case 1 : in -> primary if(modeTrack == 1) { // return track number of primary particle // of the given track kine = HGeantKine::getPrimary(numTrack,fGeantKineCat); return kine->getTrack(); } //-------------------------------------- // case 2 and 3 : entering SHOWER/TOF Int_t tempTrack = numTrack; first = 0; //-------------------------------------- // case 2 : entering SHOWER if(modeTrack == 2) { while((kine = kine->getParent(tempTrack,fGeantKineCat)) != 0) { first = kine->getFirstShowerHit(); if(first != -1) { // track is still in SHOWER HGeantShower* gshower = (HGeantShower*)fGeantShowerCat->getObject(first); Int_t s1 = gshower->getSector(); if(s == s1) { // check only sector tempTrack = kine->getTrack(); } else { // track has no SHOWER hit any longer // which fulfills the condition break; } } else { // track has no SHOWER hit any longer, // so the previous object was the one we // searched for break; } } return tempTrack; } //-------------------------------------- //-------------------------------------- // case 3 : entering TOFINO if(modeTrack >= 3) { Bool_t foundTof = kFALSE; Int_t tempTrack2 = tempTrack; do { first = kine->getFirstTofHit(); if(first != -1) { // we are in TOF HGeantTof* gtof = (HGeantTof*)fGeantTofCat->getObject(first); Int_t s1 = gtof->getSector(); Int_t m = gtof->getModule(); if(s == s1 && m > 21 ) { // check only sector + TOFINO foundTof = kTRUE; tempTrack = tempTrack2; } } tempTrack2 = kine->getParentTrack(); } while( tempTrack2 > 0 && (kine = (HGeantKine*)fGeantKineCat->getObject(tempTrack2 - 1)) != 0); if(foundTof) { tempTrack += 10000; } else { kine = (HGeantKine*)fGeantKineCat->getObject(numTrack - 1); if( modeTrack == 3 ){ // store primaries if no TOFino was found kine = HGeantKine::getPrimary(numTrack,fGeantKineCat); tempTrack = kine->getTrack() + 100000000; } else if (modeTrack == 4){ //-------------------------------------- // recover first particle entering SHOWER while((kine = kine->getParent(tempTrack,fGeantKineCat)) != 0) { first = kine->getFirstShowerHit(); if(first != -1) { // track is still in SHOWER HGeantShower* gshower = (HGeantShower*)fGeantShowerCat->getObject(first); Int_t s1 = gshower->getSector(); if(s == s1) { // check only sector tempTrack = kine->getTrack(); } else { // track has no SHOWER hit any longer // which fulfills the condition break; } } else { // track has no SHOWER hit any longer, // so the previous object was the one we // searched for break; } } tempTrack += 100000000; //-------------------------------------- } } } //-------------------------------------- return tempTrack; } Int_t HShowerPadDigitizer::findFirstHitInRpc(Int_t trackID) { //------------------------------------------------- // Used to suppress the secondaries created in the // SHOWER itself. // 0 = realistic (secondaries included) // 1 primary particle is stored // 2 the first track number entering the SHOWER in SAME SECTOR is stored // 3 the first track number entering the RPC in SAME SECTOR is stored // or the primary track if no RPC was found // 4 (default) the first track number entering the RPC in SAME SECTOR is stored // or the first track in SHOWER if no RPC was found Int_t numTrack = trackID; if(numTrack <= 0) return numTrack; // nothing to do for negative track numbers //-------------------------------------- // case 0 : in = out if(modeTrack == 0) { return numTrack; } HGeantShower *poldShower; Int_t first = 0; Int_t parent= 0; HGeantKine* kine = (HGeantKine*)fGeantKineCat->getObject(numTrack - 1); if(kine){ parent = kine->getParentTrack(); if( parent == 0) { return numTrack; } // nothing todo first = kine->getFirstShowerHit(); if(first != -1){ poldShower = (HGeantShower*)fGeantShowerCat->getObject(first); } else { Error("findFirstHitInRpc()","No first shower hit!"); return numTrack; } } else { Error("findFirstHitInRpc()","Received Zero pointer for kine!"); return numTrack; } if(numTrack != poldShower->getTrack()){ Error("findFirstHitInRpc()","First shower hit not same trackID!"); return numTrack; } //-------------------------------------------------------- // return the track number for // the selected option modeTrack Int_t s = poldShower->getSector(); //-------------------------------------- // case 1 : in -> primary if(modeTrack == 1) { // return track number of primary particle // of the given track kine = HGeantKine::getPrimary(numTrack,fGeantKineCat); return kine->getTrack(); } //-------------------------------------- // case 2 and 3 : entering SHOWER/RPC Int_t tempTrack = numTrack; first = 0; //-------------------------------------- // case 2 : entering SHOWER if(modeTrack == 2) { while((kine = kine->getParent(tempTrack,fGeantKineCat)) != 0) { first = kine->getFirstShowerHit(); if(first != -1) { // track is still in SHOWER HGeantShower* gshower = (HGeantShower*)fGeantShowerCat->getObject(first); Int_t s1 = gshower->getSector(); if(s == s1) { // check only sector tempTrack = kine->getTrack(); } else { // track has no SHOWER hit any longer // which fulfills the condition break; } } else { // track has no SHOWER hit any longer, // so the previous object was the one we // searched for break; } } return tempTrack; } //-------------------------------------- //-------------------------------------- // case 3 : entering RPC if(modeTrack >= 3) { Bool_t foundRPC = kFALSE; Int_t tempTrack2 = tempTrack; do { first = kine->getFirstRpcHit(); if(first != -1) { // we are in TOF HGeantRpc* grpc = (HGeantRpc*)fGeantRpcCat->getObject(first); Int_t s1 = grpc->getSector(); if(s == s1) { // check only sector + RPC foundRPC = kTRUE; tempTrack = tempTrack2; } } tempTrack2 = kine->getParentTrack(); } while( tempTrack2 > 0 && (kine = (HGeantKine*)fGeantKineCat->getObject(tempTrack2 - 1)) != 0); if(foundRPC) { tempTrack += 10000; } else { kine = (HGeantKine*)fGeantKineCat->getObject(numTrack - 1); if( modeTrack == 3 ){ // store primaries if no RPC was found kine = HGeantKine::getPrimary(numTrack,fGeantKineCat); tempTrack = kine->getTrack() + 100000000; } else if (modeTrack == 4){ //-------------------------------------- // recover first particle entering SHOWER while((kine = kine->getParent(tempTrack,fGeantKineCat)) != 0) { first = kine->getFirstShowerHit(); if(first != -1) { // track is still in SHOWER HGeantShower* gshower = (HGeantShower*)fGeantShowerCat->getObject(first); Int_t s1 = gshower->getSector(); if(s == s1) { // check only sector tempTrack = kine->getTrack(); } else { // track has no SHOWER hit any longer // which fulfills the condition break; } } else { // track has no SHOWER hit any longer, // so the previous object was the one we // searched for break; } } tempTrack += 100000000; //-------------------------------------- } } } //-------------------------------------- return tempTrack; } Int_t HShowerPadDigitizer::findFirstHit(Int_t trackID) { if(isRpc) {return findFirstHitInRpc (trackID);} else {return findFirstHitInShower(trackID);} } void HShowerPadDigitizer::updatePad(HShowerPad *pPad, Float_t fIndQ, Int_t nSect, Int_t nMod, Int_t nTrack) { Int_t nRow, nCol; HShowerRaw *pShowerRaw=NULL; HLocation loc; Float_t fUpdateThreshold = pDigiPar->getUpdatePadThreshold(); pPad->getPadPos(&nRow, &nCol); loc.set(4, nSect, nMod, nRow, nCol); pShowerRaw = (HShowerRaw*)((HMatrixCategory*)getOutCat())->getObject(loc); if (!pShowerRaw) { pShowerRaw = (HShowerRaw*)((HMatrixCategory*)getOutCat())->getSlot(loc); if (pShowerRaw != NULL) { pShowerRaw = new(pShowerRaw) HShowerRaw; pShowerRaw->setSector(loc[0]); pShowerRaw->setModule(loc[1]); pShowerRaw->setRow(loc[2]); pShowerRaw->setCol(loc[3]); } } if(pShowerRaw) { pShowerRaw->addCharge(fIndQ); if(pShowerRaw->getCharge() >= fUpdateThreshold) { Int_t tempTrack = findFirstHit(nTrack); //------------------------------------------------------------ Int_t addTr = pShowerRaw->getAddress(); map >::iterator pos = trackMap.find(addTr); if( pos == trackMap.end() ) { // address not yet registered vector v; v.push_back(tempTrack); trackMap.insert( make_pair( addTr, v ) ); } else { // add track to the list of the pad (pos->second).push_back(tempTrack); } //------------------------------------------------------------ } } } Float_t HShowerPadDigitizer::calcCharge(Float_t charge, Float_t dist, Float_t Xd, Float_t Yd, Float_t Xu, Float_t Yu) { //charge density matrix is calculated const Float_t twoPI = 6.28318530718; return ((charge/twoPI)*( atan(Xd*Yd/(dist*sqrt(dist*dist+Xd*Xd+Yd*Yd))) - atan(Xd*Yu/(dist*sqrt(dist*dist+Xd*Xd+Yu*Yu))) + atan(Xu*Yu/(dist*sqrt(dist*dist+Xu*Xu+Yu*Yu))) - atan(Xu*Yd/(dist*sqrt(dist*dist+Xu*Xu+Yd*Yd))))); } /* ************************************************************************** */ void HShowerPadDigitizer::digiPads(HShowerGeantWire* pWireHit) { Float_t fX, fY, fQ; Float_t factor; Float_t sum; Int_t nPadX, nPadY, nPadRange; Float_t pChargeTable[MAX_PADS_DIST]; Int_t iMaxInArray, m; HShowerPad *pTmpPad; HShowerPad *pCPad; Float_t fPadThreshold = pDigiPar->getPadThreshold(); Float_t fUpdateThreshold = pDigiPar->getUpdatePadThreshold(); HShowerPadTab *pPadParam = pGeometry->getPadParam(pWireHit->getModule()); // don't calculate if charge belowe threshold if((fQ = pWireHit->getCharge()) < fUpdateThreshold) return; pWireHit->getXY(&fX, &fY); //row = Y,col = X if((pCPad = pPadParam->getPad(fX, fY)) == NULL) return; pCPad->getPadPos(&nPadY, &nPadX); //row = Y; col = X iMaxInArray = MAX_PADS_DIST; nPadRange = (MAX_PADS_DIST - 1) / 2; if(numericalCalc(pCPad, pWireHit, pChargeTable) != 0) return; // calculate sum for the first renormalization sum = 0.0f; for(m = 0; m < iMaxInArray; m++) sum += pChargeTable[m]; if(sum <= 0.0f) return; factor = fQ / sum; for(m = 0; m < iMaxInArray; m++) { pChargeTable[m] *= factor; } // calculate sum for a renormalization sum = 0.0f; for(m = 0; m < iMaxInArray; m++) { if(pChargeTable[m] >= fPadThreshold) sum += pChargeTable[m]; else pChargeTable[m] = 0.0f; } // if sum is below threshold put all charge to the central pad if(sum < fPadThreshold) { for(m = 0; m < iMaxInArray; m++) pChargeTable[m] = 0.0f; pChargeTable[nPadRange] = fQ; } else { // renormalization: sum of the charge on the pads // must be equal to the charge on the wire factor = fQ / sum; for(m = 0; m < iMaxInArray; m++) { if(pChargeTable[m] > 0.0f) pChargeTable[m] *= factor; } } for(m = 0; m < iMaxInArray; m++) { pTmpPad = pPadParam->getPad(nPadY, nPadX + m - nPadRange); if(((m == nPadRange) || (pChargeTable[m] >= fPadThreshold)) && (pTmpPad != NULL) && (pTmpPad->getPadFlag())) { updatePad(pTmpPad, pChargeTable[m], pWireHit->getSector(), pWireHit->getModule(), pWireHit->getTrack()); } } } //--------------------------------------------------------- Double_t HShowerPadDigitizer::gainCharge(HShowerRawMatr *pRawMatr) { // gain charge in simulation return pRawMatr->getCharge()*pDigiPar->getGain(pRawMatr->getSector(), pRawMatr->getModule()); } //--------------------------------------------------------- Bool_t HShowerPadDigitizer::checkEfficiency(HShowerRawMatr *pRawMatr) { // check efficiency // model good for post converters // if pre-converter if(pRawMatr->getModule() == 0) return kTRUE; Double_t threshold = pDigiPar->getThreshold(pRawMatr->getSector(),pRawMatr->getModule()); if(pRawMatr->getCharge() > threshold) return kTRUE; else return kFALSE; } /* ************************************************************************** */ Int_t HShowerPadDigitizer::calcLimit(Float_t fCor, Int_t nMatrixRange, Float_t fBoxSize) { //charge density matrix is delimited to the sensitive area which covers fired pads Int_t k, nLimit; k = (Int_t)((fabs(fCor) + 0.5*fBoxSize)/fBoxSize); if (fCor < 0 ) nLimit = (k>nMatrixRange) ? 0 : nMatrixRange - k; else nLimit = (k>nMatrixRange) ? 2*nMatrixRange : nMatrixRange + k; return nLimit; } void HShowerPadDigitizer::moveCoord(HShowerPad *pPad,Float_t distWire, Float_t fDx, Float_t fDy, Float_t *corXld, Float_t *corYld, Float_t *corXlu, Float_t *corYlu, Float_t *corXrd, Float_t *corYrd, Float_t *corXru, Float_t *corYru) { // moves charge density matrix to the place where is placed induced charge on the sense wire *corXld = pPad->getXld() - fDx; // *corYld = pPad->getYld() - fDy; *corYld = fDy - fDy - distWire; *corXlu = pPad->getXlu() - fDx; // *corYlu = pPad->getYlu() - fDy; *corYlu = fDy - fDy + distWire; *corXrd = pPad->getXrd() - fDx; // *corYrd = pPad->getYrd() - fDy; *corYrd = fDy - fDy - distWire; *corXru = pPad->getXru() - fDx; // *corYru = pPad->getYru() - fDy; *corYru = fDy - fDy + distWire; } void HShowerPadDigitizer::analyticCalc(HShowerPad *pPad, HShowerGeantWire* pWireHit) { // this function calculates charge induced on pad analytically // it is useful only for square or rectangular pads Float_t fX, fY; Float_t fIndQ; Float_t corXld, corYld, corXlu, corYlu; Float_t corXrd, corYrd, corXru, corYru; pWireHit->getXY(&fX, &fY); Int_t nModule = pWireHit->getModule(); HShowerWireTab *pfWire = pGeometry->getWireTab(nModule); Float_t distWire = pfWire->getDistWire(); // cout << " distWire " << distWire << endl; moveCoord(pPad, distWire, fX, fY, &corXld, &corYld, &corXlu, &corYlu, &corXrd, &corYrd, &corXru, &corYru); fIndQ = calcCharge(pWireHit->getCharge(), pDigiPar->getPlaneDist(), corXld, corYld, corXru, corYru); updatePad(pPad, fIndQ, pWireHit->getSector(), pWireHit->getModule(), pWireHit->getTrack()); } //------------------------------------------------------------------------------ Int_t HShowerPadDigitizer::numericalCalc(HShowerPad *pPad, HShowerGeantWire* pWireHit, Float_t pfValues[MAX_PADS_DIST]) { Float_t fAl, fAr, fBl, fBr; Float_t fXd, fXu, fYd, fYu; Float_t fXmin, fXmax; Float_t fYmin, fYmax; Float_t fWireX, fWireY; Float_t fY, fX, fQ; Int_t i; HShowerFrameCorner *pCorner; if((pPad == NULL) || (pPad->getPadFlag() == 0)) return -1; Int_t iModule = pWireHit->getModule(); HShowerFrame *pFrame = pGeometry->getFrame(iModule); Float_t fDistWire = pGeometry->getWireTab(iModule)->getDistWire(); Int_t iMatrixRange = pDigiPar->getMatrixRange(); Float_t fBoxSize = pDigiPar->getBoxSize(); const Float_t* pfChargeMatrix = pDigiPar->getChargeMatrix(); // left border if((pPad->getPadFlag() == 2) && (pPad->getXld() < 0.0)) { pCorner = pFrame->getCorner(0); fXd = pCorner->getX(); fYd = pCorner->getY(); pCorner = pFrame->getCorner(1); fXu = pCorner->getX(); fYu = pCorner->getY(); } else { fXd = pPad->getXld(); fYd = pPad->getYld(); fXu = pPad->getXlu(); fYu = pPad->getYlu(); } fAl = (fXu - fXd) / (fYu - fYd); fBl = fXu - fAl * fYu; // right border if((pPad->getPadFlag() == 2) && (pPad->getXld() >= 0.0)) { pCorner = pFrame->getCorner(3); fXd = pCorner->getX(); fYd = pCorner->getY(); pCorner = pFrame->getCorner(2); fXu = pCorner->getX(); fYu = pCorner->getY(); } else { fXd = pPad->getXrd(); fYd = pPad->getYrd(); fXu = pPad->getXru(); fYu = pPad->getYru(); } fAr = (fXu - fXd) / (fYu - fYd); fBr = fXu - fAr * fYu; // range of the charge matrix pWireHit->getXY(&fWireX, &fWireY); fXmin = fWireX - iMatrixRange * fBoxSize; fXmax = fWireX + iMatrixRange * fBoxSize; fYmin = fWireY - fDistWire + 0.5 * fBoxSize; fYmax = fWireY + fDistWire - 0.5 * fBoxSize; memset(pfValues, 0, MAX_PADS_DIST * sizeof(Float_t)); for(fY = fYmin; fY <= fYmax; fY += fBoxSize) { // left border if((fX = fAl * fY + fBl) >= fXmin && (fWireX-fX) >= 0.0) { i = (int)((fX - fXmin) / fBoxSize); if(i>2*iMatrixRange||i<0) { cout << " left border fWireX-fX " << fWireX-fX << endl; cout << " fXd,fYd,fXu,fYu " << pPad->getXld() <<","<< pPad->getYld() << ","<< pPad->getXlu() << ","<< pPad->getYlu() << endl; Error("numericalCalc()","left border matrix : max %i, ind %i, fX %7.4f, fXmin %7.4f, fWireX %7.4f,fBoxSize %5.3f, PadNr %i",2*iMatrixRange,i,fX,fXmin,fWireX,fBoxSize,pPad->getPadNr()); continue; } pfValues[0] += pfChargeMatrix[i]; fQ = pfChargeMatrix[i]; } else fQ = 0.0f; // right border if((fX = fAr * fY + fBr) <= fXmax && (fWireX-fX) <= 0.0) { i = (int)((fX - fXmin) / fBoxSize); if(i>2*iMatrixRange||i<0) { cout << " right border fWireX-fX " << fWireX-fX << endl; cout << " fXd,fYd,fXu,fYu " << pPad->getXrd() <<","<< pPad->getYrd() << ","<< pPad->getXru() << ","<< pPad->getYru() << endl; Error("numericalCalc()","right border matrix : max %i, ind %i, fX %7.4f, fXmin %7.4f, fWireX %7.4f, fBoxSize %5.3f, PadNr %i",2*iMatrixRange,i,fX,fXmin,fWireX,fBoxSize,pPad->getPadNr()); continue; } pfValues[1] += pfChargeMatrix[i] - fQ; pfValues[2] += pfChargeMatrix[2 * iMatrixRange] - pfChargeMatrix[i]; } else pfValues[1] += pfChargeMatrix[2 * iMatrixRange] - fQ; } return 0; }