00001
00002 #include "StMuDSTMaker/COMMON/StMuDst.h"
00003 #include "StMuDSTMaker/COMMON/StMuEvent.h"
00004 #include "StMuDSTMaker/COMMON/StMuDstMaker.h"
00005 #include "StMuDSTMaker/COMMON/StMuEmcCollection.h"
00006 #include "StMuDSTMaker/COMMON/StMuEmcUtil.h"
00007 #include "StMuDSTMaker/COMMON/StMuTrack.h"
00008 #include "StMuDSTMaker/COMMON/StMuPrimaryVertex.h"
00009
00010 #include "StEvent/StTriggerId.h"
00011 #include "StEvent/StEvent.h"
00012 #include "StEventTypes.h"
00013
00014 #include "StDetectorDbMaker/StDetectorDbTriggerID.h"
00015 #include "StEvent/StL0Trigger.h"
00016
00017 #include "StDaqLib/TRG/trgStructures2005.h"
00018
00019
00020 #include "StEmcUtil/database/StBemcTables.h"
00021 #include "StEmcUtil/database/StEmcDecoder.h"
00022 #include "StEmcUtil/geometry/StEmcGeom.h"
00023 #include "StEmcUtil/projection/StEmcPosition.h"
00024 #include "StEmcADCtoEMaker/StBemcData.h"
00025 #include "StEmcADCtoEMaker/StEmcADCtoEMaker.h"
00026 #include "StEmcRawMaker/StBemcRaw.h"
00027 #include "StEmcRawMaker/defines.h"
00028 #include "StEmcRawMaker/StBemcTables.h"
00029 #include "StEmcTriggerMaker/StEmcTriggerMaker.h"
00030
00031
00032 #include "StMessMgr.h"
00033
00034
00035 #include "StTriggerStudyMaker.h"
00036 #include "StDetectorDbMaker/StDetectorDbBeamInfo.h"
00037 #include "StSpinPool/StTriggerFilterMaker/StTriggerFilterMaker.h"
00038 #include "StTriggerStudyEvent.h"
00039 #include "StEnumerations.h"
00040 #include "StTriggerUtilities/StTriggerSimuMaker.h"
00041 #include "StTriggerUtilities/StTriggerSimuResult.h"
00042 #include "TTree.h"
00043 #include "TFile.h"
00044
00045 ClassImp(StTriggerStudyMaker)
00046
00047 StTriggerStudyMaker::StTriggerStudyMaker(const char* filename, const char* name):StMaker(name),nTriggers(12)
00048 {
00049 mName = filename;
00050 runNumber = 0;
00051 muDstMaker = NULL;
00052 mTriggerSimuMaker = NULL;
00053 mTree = NULL;
00054 mFile = NULL;
00055 mTSEvent = NULL;
00056 }
00057
00058 Int_t StTriggerStudyMaker::Init()
00059 {
00060 int trigarray[12] = {117001,137221,137222,137822,117300,137571,137575,137585,137586,137611,137501,137622};
00061 for(int i = 0; i < nTriggers; i++){
00062 triggers[i] = trigarray[i];
00063 trigmap[triggers[i]] = i;
00064 }
00065
00066 mFile = new TFile(mName,"RECREATE");
00067 mTSEvent = new StTriggerStudyEvent();
00068 mTree = new TTree("trigStudyTree","Trigger Study Tree");
00069 mTree->Branch("event_branch","StTriggerStudyEvent",&mTSEvent);
00070 mTree->SetAutoSave(1000000000);
00071
00072 muDstMaker = dynamic_cast<StMuDstMaker*>(GetMaker("MuDst"));
00073 assert(muDstMaker);
00074 mADCtoEMaker = dynamic_cast<StEmcADCtoEMaker*>(GetMaker("Eread")); assert(mADCtoEMaker);
00075 mTriggerSimuMaker = dynamic_cast<StTriggerSimuMaker*>(GetMaker("StarTrigSimu"));
00076
00077 mDecoder = new StEmcDecoder();
00078 mEmcGeom = StEmcGeom::instance("bemc");
00079 mBemcTables = mADCtoEMaker->getBemcData()->getTables();
00080
00081 zvertall = new TH1F("zvertall","Z Vertex Distribution",280,-139.5,140.5);
00082 bbctall = new TH1F("bbctall","BBC Time Bins",18,-1.5,16.5);
00083 jp1et = new TH1F("jp1et","JP1 Et Distribution",70,17.9,31.9);
00084 jp2et = new TH1F("jp2et","JP2 Et Distribution",70,17.9,31.9);
00085 jp1et0 = new TH1F("jp1et0","JP1 Et Distribution",70,3.9,17.9);
00086 jp0et0 = new TH1F("jp0et0","JP0 Et Distribution",70,3.9,17.9);
00087 for(int i = 0; i < nTriggers; i++){
00088 trigmap[triggers[i]] = i;
00089 char namez[100];
00090 sprintf(namez,"zvert%i",triggers[i]);
00091 zverttrig[i] = new TH1F(namez,"Z Vertex Distribution",280,-139.5,140.5);
00092 char nameb[100];
00093 sprintf(nameb,"bbct%i",triggers[i]);
00094 bbcttrig[i] = new TH1F(nameb,"BBC Time Bins",18,-1.5,16.5);
00095 }
00096
00097 return StMaker::Init();
00098 }
00099
00100 Int_t StTriggerStudyMaker::InitRun(int run)
00101 {
00102 if(runNumber != 0){
00103 LOG_ERROR<<"Attempting to change run: check file list"<<endm;
00104 return kStErr;
00105 }
00106 LOG_INFO<<"Now starting Run#"<<run<<endm;
00107 runNumber = run;
00108 return kStOk;
00109
00110 }
00111
00112 Int_t StTriggerStudyMaker::Finish(){
00113 mFile->cd();
00114 mTree->Write();
00115 zvertall->Write();
00116 bbctall->Write();
00117 jp1et->Write();
00118 jp2et->Write();
00119 jp1et0->Write();
00120 jp0et0->Write();
00121 for(int i = 0; i < nTriggers; i++){
00122 zverttrig[i]->Write();
00123 bbcttrig[i]->Write();
00124 }
00125 mFile->Close();
00126 return StMaker::Finish();
00127 }
00128
00129 Int_t StTriggerStudyMaker::Make(){
00130 fillTree();
00131 mTSEvent->Clear();
00132 return kStOK;
00133 }
00134
00135 void StTriggerStudyMaker::fillTree()
00136 {
00137 StDetectorDbTriggerID* dbtrig = StDetectorDbTriggerID::instance();
00138 StMuDst* mudst = muDstMaker->muDst();assert(mudst);
00139 StMuEvent* muevent = mudst->event();assert(muevent);
00140
00141 int bbctimebin = muevent->bbcTriggerDetector().onlineTimeDifference();
00142 mTSEvent->setBbcTimeBin(bbctimebin);
00143 bbctall->Fill(bbctimebin/32);
00144
00145 float vertmag = TMath::Abs(muevent->primaryVertexPosition().mag());
00146 float vertz = muevent->primaryVertexPosition().z();
00147 if(vertmag > 1e-7)mTSEvent->setVertexPosition(muevent->primaryVertexPosition());
00148 if(vertmag > 1e-7)zvertall->Fill(vertz);
00149
00150 int isjp1 = 0;
00151 int isjp0 = 0;
00152 int isjp2 = 0;
00153 float jpEt[12] = {0,0,0,0,0,0,0,0,0};
00154 float maxjpEt = -1;
00155 mEmcCollection = mudst->emcCollection();
00156 StDetectorId detectorid = static_cast<StDetectorId>(kBarrelEmcTowerId);
00157 StEmcDetector* detector=mEmcCollection->detector(detectorid);
00158 if(detector){
00159 for(int m=1;m<=120;m++){
00160 StEmcModule* module = detector->module(m);
00161 if(module){
00162 StSPtrVecEmcRawHit& rawHit=module->hits();
00163 for(unsigned int k=0;k<rawHit.size();k++){
00164 if(rawHit[k]){
00165 int module = rawHit[k]->module();
00166 int eta = rawHit[k]->eta();
00167 int submodule = abs(rawHit[k]->sub());
00168 int ADC = rawHit[k]->adc();
00169 int hitid;
00170 int stat = mEmcGeom->getId(module,eta,submodule,hitid);
00171 int status;
00172 mBemcTables->getStatus(BTOW,hitid,status);
00173 float pedestal,rms;
00174 mBemcTables->getPedestal(BTOW,hitid,0,pedestal,rms);
00175 if(stat==0&&status==1&&((float)ADC-pedestal) > 2*rms){
00176
00177 float theta;
00178 mEmcGeom->getTheta(hitid,theta);
00179 float calib = mBemcTables->calib(BTOW,hitid);
00180 float energy = calib * (float)(ADC-pedestal);
00181 float et = energy*sin(theta);
00182 if(et > 0.2){
00183 int jpID;
00184 mDecoder->GetJetPatchFromTowerId(hitid,jpID);
00185 jpEt[jpID]+=et;
00186 }
00187 }
00188 }
00189 }
00190 }
00191 else { LOG_WARN<<"couldn't find StEmcModule "<<m<<" for detector "<<BTOW<<endm;}
00192 }
00193 }
00194
00195 const StTriggerId& trigs = muevent->triggerIdCollection().nominal();
00196 vector<unsigned int>triggers = trigs.triggerIds();
00197 vector<unsigned int>prescales;
00198 vector<unsigned int>simutriggers;
00199
00200 mTSEvent->setTriggers(triggers);
00201
00202 for(unsigned int i = 0; i < triggers.size(); i++){
00203 unsigned int prs = (unsigned int)dbtrig->getTotalPrescaleByTrgId(triggers[i]);
00204 prescales.push_back(prs);
00205
00206 if(triggers[i] == 137222)isjp1 = 1;
00207 if(triggers[i] == 137501)isjp0 = 1;
00208 if(triggers[i] == 137585)isjp2 = 1;
00209 map<int,int>::iterator iter = trigmap.find(triggers[i]);
00210 if(iter!=trigmap.end()){
00211 int tr = iter->second;
00212 if(vertmag > 1e-7)zverttrig[tr]->Fill(vertz);
00213 bbcttrig[tr]->Fill(bbctimebin/32);
00214 }
00215 if(mTriggerSimuMaker && mTriggerSimuMaker->isTrigger(triggers[i])){
00216 simutriggers.push_back(triggers[i]);
00217 StTriggerSimuResult trigsimuresult = mTriggerSimuMaker->detailedResult(triggers[i]);
00218 vector<short> towers = trigsimuresult.highTowerIds();
00219 for(unsigned int j = 0; j < towers.size(); j++){
00220 int dsm = trigsimuresult.highTowerAdc(towers[j]);
00221
00222 }
00223 vector<short> jpatches = trigsimuresult.jetPatchIds();
00224 for(unsigned int j = 0; j < jpatches.size(); j++){
00225 int dsm = trigsimuresult.jetPatchAdc(jpatches[j]);
00226 if(dsm > 0)mTSEvent->addJPatchDSM(jpatches[j],dsm);
00227 }
00228 vector<short> tpatches = trigsimuresult.triggerPatchIds();
00229 for(unsigned int j = 0; j < tpatches.size(); j++){
00230 int dsm = trigsimuresult.triggerPatchAdc(tpatches[j]);
00231 if(dsm > 0)mTSEvent->addTPatchDSM(tpatches[j],dsm);
00232 }
00233 }
00234 }
00235 vector<float> jpetvec;
00236 for(int k = 0; k < 12; k++){
00237 if(jpEt[k] > maxjpEt)maxjpEt = jpEt[k];
00238 jpetvec.push_back(jpEt[k]);
00239 }
00240 if(isjp1) jp1et->Fill(maxjpEt);
00241 if(isjp1 && isjp2)jp2et->Fill(maxjpEt);
00242 if(isjp0) jp0et0->Fill(maxjpEt);
00243 if(isjp0 && isjp1)jp1et0->Fill(maxjpEt);
00244
00245 mTSEvent->setJPEt(jpetvec);
00246 mTSEvent->setPrescales(prescales);
00247 mTSEvent->setSimuTriggers(simutriggers);
00248
00249 mTree->Fill();
00250 }
