00001 #include "StMessMgr.h"
00002 #include "StLuminosityMaker.h"
00003 #include "StLuminosityHolder.h"
00004 #include "StDetectorDbMaker/StDetectorDbTriggerID.h"
00005 #include "StDetectorDbMaker/StDetectorDbBeamInfo.h"
00006 #include "StMuDSTMaker/COMMON/StMuDst.h"
00007 #include "StMuDSTMaker/COMMON/StMuEvent.h"
00008 #include "StEvent/StTriggerId.h"
00009 #include "StSpinPool/StTriggerFilterMaker/StTriggerFilterMaker.h"
00010
00011 ClassImp(StLuminosityMaker)
00012
00013 StLuminosityMaker::StLuminosityMaker(const char* name):StMaker(name)
00014 {
00015 mTriggers.clear();
00016 mNTotal.clear();
00017 mNCuts.clear();
00018 mPrescales.clear();
00019 mLumTotal.clear();
00020 mLumCuts.clear();
00021 mXsec = -1.0;
00022 runMode = "";
00023 runNumber = 0;
00024 mVertexCut = 0;
00025 mOverrideMode = 0;
00026 muDstMaker = NULL;
00027 mTriggerSimuMaker = NULL;
00028 mLumHolder = NULL;
00029 mFilterMode = 0;
00030 }
00031
00032 Int_t StLuminosityMaker::Init()
00033 {
00034 muDstMaker = dynamic_cast<StMuDstMaker*>(GetMaker("MuDst"));
00035 assert(muDstMaker);
00036 mTriggerSimuMaker = dynamic_cast<StTriggerSimuMaker*>(GetMaker("StarTrigSimu"));
00037 if(!strcmp(runMode,"")){
00038 LOG_ERROR<<"runMode is not set"<<endm;
00039 return kStErr;
00040 }
00041
00042 LOG_INFO<<"Using minbias cross-section: "<<getCrossSectionNB()<<" (nb)"<<endm;
00043
00044 mLumHolder = new TClonesArray("StLuminosityHolder",20);
00045
00046 return StMaker::Init();
00047 }
00048
00049 Int_t StLuminosityMaker::InitRun(int run)
00050 {
00051 if(runNumber != 0){
00052 saveOutput();
00053 printOutput();
00054 }
00055 LOG_INFO<<"Now starting Run#"<<run<<endm;
00056 runNumber = run;
00057 for(unsigned int i = 0; i < mTriggers.size(); i++){
00058 mNTotal[i] = 0;
00059 mNCuts[i] = 0;
00060 mNVertex[i] = 0;
00061 mPrescales[i] = 0.0;
00062 }
00063
00064 return kStOk;
00065
00066 }
00067
00068 Int_t StLuminosityMaker::Finish(){
00069 saveOutput();
00070 printOutput();
00071 return StMaker::Finish();
00072 }
00073
00074 Int_t StLuminosityMaker::Make(){
00075 StDetectorDbTriggerID* dbtrig = StDetectorDbTriggerID::instance();
00076 StMuDst* mudst = muDstMaker->muDst();assert(mudst);
00077 StMuEvent* muevent = mudst->event();assert(muevent);
00078
00079 int passtrig = 0;
00080 int passvert = 0;
00081
00082 const StTriggerId& trigs = muevent->triggerIdCollection().nominal();
00083
00084 for(unsigned int i = 0; i < mTriggers.size(); i++){
00085 if(trigs.isTrigger(mTriggers[i])){
00086 float prs = dbtrig->getTotalPrescaleByTrgId(mTriggers[i]);
00087 mPrescales[i] = prs;
00088 mNTotal[i]++;
00089 passtrig++;
00090 float vertz = TMath::Abs(muevent->primaryVertexPosition().z());
00091 if(vertz > 1e-5)mNVertex[i]++;
00092 if(mTriggerSimuMaker && !mTriggerSimuMaker->isTrigger(mTriggers[i]))continue;
00093 mNSoftTrig[i]++;
00094 if(mVertexCut > vertz && vertz > 1e-5){
00095 mNCuts[i]++;
00096 passvert++;
00097 }
00098 }
00099 }
00100 if(mFilterMode && !passvert && !passtrig) return kStSkip;
00101
00102 return kStOK;
00103 }
00104
00105 void StLuminosityMaker::addTrigger(unsigned int trigId)
00106 {
00107 mTriggers.push_back(trigId);
00108 mNTotal.push_back(0);
00109 mNCuts.push_back(0);
00110 mNVertex.push_back(0);
00111 mNSoftTrig.push_back(0);
00112 mPrescales.push_back(0.0);
00113 mLumTotal.push_back(0.0);
00114 mLumCuts.push_back(0.0);
00115 mLumVertex.push_back(0.0);
00116 mLumSoftTrig.push_back(0.0);
00117 }
00118
00119 void StLuminosityMaker::setMode(char* newMode)
00120 {
00121 float oldxsec = mXsec;
00122
00123 if(!strcmp(newMode,"dAu2008")){
00124 runMode = newMode;
00125 mXsec = 13.8e3;
00126 }else if(!strcmp(newMode,"pp2008")){
00127 runMode = newMode;
00128 mXsec = 26.1e6;
00129 }else if(!strcmp(newMode,"AuAu2007")){
00130 runMode = newMode;
00131 mXsec = 10e9;
00132 }else if(!strcmp(newMode,"pp2006")){
00133 runMode = newMode;
00134 mXsec = 26.1e6;
00135 }else if(!strcmp(newMode,"pp2006_62GeV")){
00136 runMode = newMode;
00137 mXsec = 20e6;
00138 }else if(!strcmp(newMode,"pp2005")){
00139 runMode = newMode;
00140 mXsec = 26.1e6;
00141 }else if(!strcmp(newMode,"CuCu2005")){
00142 runMode = newMode;
00143 mXsec = 3.3e9;
00144 }else if(!strcmp(newMode,"pp2004")){
00145 runMode = newMode;
00146 mXsec = 26.1e6;
00147 }else if(!strcmp(newMode,"AuAu2004")){
00148 runMode = newMode;
00149 mXsec = 6.5e9;
00150 }else if(!strcmp(newMode,"dAu2003")){
00151 runMode = newMode;
00152 mXsec = 2.2e9 * .2;
00153 }else if(!strcmp(newMode,"pp2003")){
00154 runMode = newMode;
00155 mXsec = 26.1e6;
00156 }else if(strcmp(newMode,"")){
00157 LOG_ERROR<<"Invalid mode "<<newMode<<" passed to StLuminosityMaker"<<endm;
00158 }else{
00159 LOG_ERROR<<"StLuminosityMaker requires a running mode"<<endm;
00160 }
00161 if(mOverrideMode){
00162 mXsec = oldxsec;
00163 }
00164
00165 }
00166
00167 void StLuminosityMaker::setCrossSectionNB(float newXsec,int overRideMode)
00168 {
00169 mXsec = newXsec;
00170 mOverrideMode = overRideMode;
00171 runMode = "manual";
00172 }
00173
00174 void StLuminosityMaker::printOutput()
00175 {
00176 LOG_INFO<<"Luminosity analyzed for Run #"<<runNumber<<":"<<endm;
00177 for(unsigned int i = 0; i < mTriggers.size(); i++){
00178 float totalsampled = 0;
00179 float cutsampled = 0;
00180 if(mNTotal[i] > 0){
00181 totalsampled = mPrescales[0]*mNTotal[0]/(mPrescales[i]*mXsec);
00182 cutsampled = totalsampled * mNCuts[i] / mNTotal[i];
00183 }
00184
00185 LOG_INFO<<"Trigger "<<mTriggers[i]<<": "<<totalsampled<<" nb^-1 total and "<<cutsampled<<" nb^-1 after cuts, "<<mNTotal[i]<<" "<<mNCuts[i]<<endm;
00186 }
00187 }
00188
00189 void StLuminosityMaker::getTriggersFromFilterMaker(const char* filtername)
00190 {
00191 StTriggerFilterMaker* filter = dynamic_cast<StTriggerFilterMaker*>(GetMaker(filtername));
00192 assert(filter);
00193 vector<unsigned int>filtertriggers = filter->getTriggers();
00194 for(unsigned int i = 0; i < filtertriggers.size(); i++){
00195 addTrigger(filtertriggers[i]);
00196 }
00197 }
00198
00199 void StLuminosityMaker::saveOutput()
00200 {
00201 for(unsigned int i = 0; i < mTriggers.size(); i++){
00202 float totalsampled = 0;
00203 float vertsampled = 0;
00204 float softtrigsampled = 0;
00205 float cutsampled = 0;
00206 if(mNTotal[i] > 0){
00207 totalsampled = mPrescales[0]*mNTotal[0]/(mPrescales[i]*mXsec);
00208 cutsampled = totalsampled * mNCuts[i] / mNTotal[i];
00209 softtrigsampled = totalsampled * mNSoftTrig[i] / mNTotal[i];
00210 vertsampled = totalsampled * mNVertex[i] / mNTotal[i];
00211 }
00212 mLumTotal[i] = totalsampled;
00213 mLumCuts[i] = cutsampled;
00214 mLumSoftTrig[i] = softtrigsampled;
00215 mLumVertex[i] = vertsampled;
00216 }
00217
00218 StLuminosityHolder* lum = new StLuminosityHolder(runNumber);
00219
00220 lum->setTriggers(mTriggers);
00221 lum->setLumTotal(mLumTotal);
00222 lum->setLumCuts(mLumCuts);
00223 lum->setLumVertex(mLumVertex);
00224 lum->setLumSoftTrig(mLumSoftTrig);
00225 lum->setNTotal(mNTotal);
00226 lum->setNCuts(mNCuts);
00227 lum->setNVertex(mNVertex);
00228 lum->setNSoftTrig(mNSoftTrig);
00229 lum->setPrescales(mPrescales);
00230 lum->setCrossSectionNB(mXsec);
00231 lum->setVertexCutcm(mVertexCut);
00232
00233 new ( (*mLumHolder)[mLumHolder->GetLast()+1] ) StLuminosityHolder(*lum);
00234
00235 }
