StRoot  1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
StSsdEmbeddingMaker.cxx
1  /***************************************************************************
2  *
3  * Author: Bouchet Jonathan
4  ***************************************************************************
5  *
6  * Description: SsdEmbeddingMaker class
7  *
8  **************************************************************************/
9 
10 #include "Stiostream.h"
11 #include <cmath>
12 using namespace std;
13 
14 #include "TH2.h"
15 #include "TFile.h"
16 #include "TNtuple.h"
17 
18 #include "St_DataSetIter.h"
19 #include "St_ObjectSet.h"
20 #include "StSsdEmbeddingMaker.h"
21 #include "StMessMgr.h"
22 
23 #include "tables/St_sls_strip_Table.h"
24 #include "tables/St_spa_strip_Table.h"
25 #include "StSsdDbMaker/StSsdDbMaker.h"
26 #include "StSsdUtil/StSsdBarrel.hh"
27 #include "StDAQMaker/StDAQReader.h"
28 #include "StDAQMaker/StSSDReader.h"
29 #include "StSsdSimulationMaker/St_spa_Maker.h"
30 #include "StSsdPointMaker/StSsdPointMaker.h"
31 
32 ClassImp(StSsdEmbeddingMaker)
33 
34 //____________________________________________________________________________
35 StSsdEmbeddingMaker::StSsdEmbeddingMaker(const char *name):StMaker(name)
36 {
38  mDoEmbedding = kTRUE;
39  mPlainSimIfNoSSD = kFALSE;
41  mSsdrealData = 0;
42  mSsdSimuData = 0;
43  }
44 //____________________________________________________________________________
45 StSsdEmbeddingMaker::~StSsdEmbeddingMaker()
46 {
47 }
48 //____________________________________________________________________________
49 Int_t StSsdEmbeddingMaker::Init()
50 {
51  if(IAttr(".histos")){
52  hStripsSimu = new TH1F("hStripsSimu","Log_{10} of the number of strip after dAQ cut in simulation",100,0,10);
53  hStripsReal = new TH1F("hStripsReal","Log_{10} of the number of strip after dAQ cut in real data",100,0,10);
54  hStripsCommon = new TH1F("hStripsCommon","number of strips in common from real file and simu file",1000,0,1000);
55  }
56  return StMaker::Init();
57 }
58 //____________________________________________________________________________
60 Int_t StSsdEmbeddingMaker::InitRun(int RunNo)
61 {
62  assert(StSsdBarrel::Instance());
63  TDataSet *ssdparams = GetInputDB("Geometry/ssd");
64  if (! ssdparams) {
65  LOG_ERROR << "No access to Geometry/ssd parameters" << endm;
66  return kStFatal;
67  }
68  TDataSetIter local(ssdparams);
69  m_dimensions = (St_ssdDimensions *)local("ssdDimensions");
70  ssdDimensions_st *dimensions = m_dimensions->GetTable();
71  setSsdParameters(dimensions);
72  return kStOk;
73 }
74 //____________________________________________________________________________
75 Int_t StSsdEmbeddingMaker::Make()
76 {
77  LOG_INFO<<"##################################"<<endm;
78  LOG_INFO<<"#### START OF SSD EMBEDDING ####"<<endm;
79  //now decide in which mode to realy run, based on options
80  //this could be done in InitRun, but could there be mixed event in one run with SSD and without SSD, ie.. some trigger mix
81  mRunningEmbedding=mDoEmbedding;
82 
83  if (mDoEmbedding && NoSsd() ){//check if SSD is present if not skip simulation
84 
85  if (mPlainSimIfNoSSD){//run plain simulation instead of embedding
86  mRunningEmbedding=kFALSE; //run plain simulation
87  }
88  else
89  { //clear data and get out
90  LOG_INFO <<"SSD SlowSimulation: SKIPPING THIS EVENT - no SSD in real data!!"<<endm;
91  return kStOk;
92  }
93  }
94  Int_t res;
95  res = GetSsdData();
96  if (res!=kStOk) return res;
97  //write out the state of simulation
98  char st[20];
99  if (mRunningEmbedding) sprintf(st,"EMBEDDING");
100  else sprintf(st,"PLAIN SIMULATION");
101  LOG_INFO<<"SSD SlowSimulation is running in the state of :"<<st<<endm;
102  if (mRunningEmbedding){ //do the embedding
103  MODE = 1;
104  Int_t SIZEREAL = mSsdrealData->GetNRows();
105  Int_t SIZESIMU = mSsdSimuData->GetNRows();
106  LOG_INFO <<"#### NUMBER OF STRIPS (SIMU) : "<<SIZESIMU <<" ####"<<endm;
107  LOG_INFO <<"#### NUMBER OF STRIPS (REAL) : "<<SIZEREAL <<" ####"<<endm;
108  if(Debug())CheckTables();
109  Int_t numberofstrips_in_common = AddRawData();
110  LOG_INFO <<"#### MERGE SIMU+REAL : "<<numberofstrips_in_common << "####" << endm;
111  if (IAttr(".histos")){
112  if(SIZESIMU>0){hStripsSimu->Fill(TMath::Log10(SIZESIMU));}
113  if(SIZEREAL>0){hStripsReal->Fill(TMath::Log10(SIZEREAL));}
114  hStripsCommon->Fill(numberofstrips_in_common);
115  }
116  }
117  return kStOk;
118 }
119 //____________________________________________________________________________
120 Int_t StSsdEmbeddingMaker::GetSsdData()
121 {
122  //EmbeddingMaker requires some data(at least empty) from the SimulationMaker
123  if(!mRunningEmbedding) return kStOk;//dont read real data if you don't need them
124  mSsdSimuData = (St_spa_strip*)GetDataSet("SsdSimuData");
125  //mSsdSimuData = (St_spa_strip*)GetDataSet("spa_strip/.data/spa_strip");
126  if (!mSsdSimuData || mSsdSimuData->GetNRows()==0){
127  LOG_ERROR << "GetSsdData : no simu input (fired strip for the SSD)"<<endm;
128  return kStErr;
129  }
130  if(mSsdSimuData) LOG_INFO << "Numbers of strips (simu) = " << mSsdSimuData->GetNRows()<< endm;
131 
132  mSsdrealData = (St_spa_strip*)GetDataSet("SsdRealData");
133  if (!mSsdrealData || mSsdrealData->GetNRows()==0){
134  LOG_ERROR << "GetSsdData : no real input (fired strip for the SSD)"<<endm;
135  return kStErr;
136  }
137  if(mSsdrealData) LOG_INFO << "Numbers of strips (real) = " << mSsdrealData->GetNRows()<< endm;
138  return kStOk;
139 }
140 //____________________________________________________
141 Int_t StSsdEmbeddingMaker::AddRawData(){
142  mSsdSimuReal = new St_spa_strip("spa_strip",50000);
143  AddData(mSsdSimuReal);
144  spa_strip_st out_strip ;
145  Int_t idWafSimu = 0;
146  Int_t iWafSimu = 0;
147  Int_t iLadSimu = 0;
148  Int_t nStripSimu = 0;
149  Int_t iSideSimu = 0;
150  spa_strip_st *strip_simu = mSsdSimuData->GetTable();
151  Int_t idWafReal = 0;
152  Int_t iWafReal = 0;
153  Int_t iLadReal = 0;
154  Int_t nStripReal = 0;
155  Int_t iSideReal = 0;
156  spa_strip_st *strip_real = mSsdrealData->GetTable();
157  Int_t currRecord = 0;
158  Int_t currRecordSimu = 0;
159  Int_t currRecordReal = 0;
160  Int_t currRecordCommon = 0;
161  Int_t IdCommon[mSsdSimuData->GetNRows()];
162  for(Int_t ll=0;ll<mSsdSimuData->GetNRows();ll++) IdCommon[ll]=0;
163  for (Int_t i = 0 ; i < mSsdrealData->GetNRows(); i++)
164  {
165  nStripReal = (int)(strip_real[i].id_strip/100000.);
166  idWafReal = strip_real[i].id_strip-10000*((int)(strip_real[i].id_strip/10000.));
167  iWafReal = idWaferToWafer(idWafReal);
168  iLadReal = (int)(idWafReal - mSsdLayer*1000 - (iWafReal+1)*100 - 1);
169  iSideReal = (strip_real[i].id_strip - nStripReal*100000 - idWafReal)/10000;
170  Int_t in =0 ;
171  for (Int_t j = 0 ; j < mSsdSimuData->GetNRows(); j++)
172  {
173  nStripSimu = (int)(strip_simu[j].id_strip/100000.);
174  idWafSimu = strip_simu[j].id_strip-10000*((int)(strip_simu[j].id_strip/10000.));
175  iWafSimu = idWaferToWafer(idWafSimu);
176  iLadSimu = (int)(idWafSimu - mSsdLayer*1000 - (iWafSimu+1)*100 - 1);
177  iSideSimu = (strip_simu[j].id_strip - nStripSimu*100000 - idWafSimu)/10000;
178 
179  if((nStripReal == nStripSimu) && (idWafReal == idWafSimu) && (iSideReal == iSideSimu)) // we found the same strip in simu and real
180  {
181  LOG_DEBUG<<Form("simu side=%d idWafer=%d Ladder=%d wafer=%d nstrip=%d signal=%d",iSideSimu,idWafSimu,iLadSimu,iWafSimu,nStripSimu,strip_simu[j].adc_count)<<endm;
182  LOG_DEBUG<<Form("real side=%d idWafer=%d Ladder=%d wafer=%d nstrip=%d signal=%d",iSideReal,idWafReal,iLadReal,iWafReal,nStripReal,strip_real[i].adc_count)<<endm;
183  out_strip.id = currRecord + 1;
184  out_strip.adc_count = (int)(strip_simu[j].adc_count + strip_real[i].adc_count);
185  out_strip.id_strip = 10000*(10*nStripSimu + iSideSimu)+idWafSimu;
186  for (Int_t e = 0 ; e < 5 ; e++)
187  {
188  out_strip.id_mchit[e] = strip_simu[j].id_mchit[e];
189  out_strip.id_mctrack[e] = strip_simu[j].id_mctrack[e];
190  }
191  mSsdSimuReal->AddAt(&out_strip);
192  IdCommon[currRecordCommon] = (int)strip_real[i].id_strip; // we keep this id
193  currRecord++;
194  currRecordCommon++;
195  in = 1;
196  LOG_DEBUG<<Form("# of strips in common =%d id=%d",currRecordCommon,(int)strip_real[i].id_strip)<<endm;
197  if(Debug())
198  {
199  for(Int_t ii=0;ii<5;ii++){
200  LOG_INFO <<"mchit["<<ii<<"]="<<strip_simu[j].id_mchit[ii]<<endm;
201  LOG_INFO <<"mctrack["<<ii<<"]="<<strip_simu[j].id_mctrack[ii]<<endm;
202  LOG_INFO << " "<< endm;
203  }
204  }
205  }
206  break; // break the loop
207  }
208  if(in==0) //no common strips, now we fill the real strip
209  {
210  out_strip.id = currRecord + 1;
211  out_strip.adc_count = (int)(strip_real[i].adc_count);
212  out_strip.id_strip = 10000*(10*nStripReal + iSideReal)+idWafReal;
213  for (Int_t e = 0 ; e < 5 ; e++)
214  {
215  out_strip.id_mchit[e] = 0;
216  out_strip.id_mctrack[e] = 0;
217  }
218  mSsdSimuReal->AddAt(&out_strip);
219  currRecordReal++;
220  currRecord++;
221  }
222  }
223  // we continue to fill the table with simu strips not found (in common with the real)
224  Int_t found ;
225  for(Int_t k=0; k <mSsdSimuData->GetNRows();k++)
226  {
227  found = 0;
228  for (Int_t kk=0;kk<currRecordCommon;kk++)
229  {
230  if((int)strip_simu[k].id_strip == IdCommon[kk]){
231  LOG_DEBUG<<Form("k = %d strip id =%d",k,(int)strip_simu[k].id_strip)<<endm;
232  found = 1;
233  }
234  }
235  if(found==0){
236  nStripSimu = (int)(strip_simu[k].id_strip/100000.);
237  idWafSimu = strip_simu[k].id_strip-10000*((int)(strip_simu[k].id_strip/10000.));
238  iWafSimu = idWaferToWafer(idWafSimu);
239  iLadSimu = (int)(idWafSimu - mSsdLayer*1000 - (iWafSimu+1)*100 - 1);
240  iSideSimu = (strip_simu[k].id_strip - nStripSimu*100000 - idWafSimu)/10000;
241  out_strip.id = currRecord + 1;
242  out_strip.adc_count = (int)(strip_simu[k].adc_count);
243  out_strip.id_strip = 10000*(10*nStripSimu + iSideSimu)+idWafSimu;
244  for (Int_t e = 0 ; e < 5 ; e++)
245  {
246  out_strip.id_mchit[e] = strip_simu[k].id_mchit[e];
247  out_strip.id_mctrack[e] = strip_simu[k].id_mctrack[e];
248  }
249  mSsdSimuReal->AddAt(&out_strip);
250  currRecord++;
251  currRecordSimu++;
252  }
253  }
254  LOG_INFO<<Form("# in common = %d totstrips after add = %d totstrips for simu only = %d for real only = %d",currRecordCommon,currRecord,currRecordSimu,currRecordReal)<<endm;
255  return currRecord;
256 }
257 //____________________________________________________________________________
258 Int_t StSsdEmbeddingMaker::Finish()
259 {
260  return kStOK;
261 }
262 //____________________________________________________________________________
263 void StSsdEmbeddingMaker::setDoEmbedding(Bool_t doIt){
264  mDoEmbedding = doIt;
265 }
266 //____________________________________________________________________________
267 void StSsdEmbeddingMaker::setPlainSimEvenIfNoSSD(Bool_t doIt){
268  mPlainSimIfNoSSD= doIt;
269 }
270 //____________________________________________________________________________
271 Int_t StSsdEmbeddingMaker::NoSsd()
272 {
273  St_DataSet *dataSet;
274 
275  dataSet = GetDataSet("StDAQReader");
276  if(!dataSet){
277  LOG_ERROR<<("BIG TROUBLE: cannot find StDAQReader in the chain and you want to run embedding")<<endm;
278  return kTRUE;
279  }
280  StDAQReader* daqReader = (StDAQReader*)(dataSet->GetObject());
281  if (!daqReader){
282  LOG_ERROR<<("BIG TROUBLE: StDAQReader is empty and you want to run embedding")<<endm;
283  return kTRUE;
284  }
285  if (!daqReader->SSDPresent ())
286  {
287  LOG_INFO<<("NO SSD in DAQ")<<endm;
288  return kTRUE; //No SSD in the datastream
289  }
290 
291  LOG_INFO<<("SSD found in DAQ")<<endm;
292  return kFALSE;
293 }
294 //______________________________________________________________________
295 void StSsdEmbeddingMaker::CheckTables(){
296  spa_strip_st *strip = mSsdSimuData->GetTable();
297 
298  Int_t idWaf = 0;
299  Int_t iWaf = 0;
300  Int_t iLad = 0;
301  Int_t nStrip = 0;
302  Int_t iSide = 0;
303  Int_t i=0;
304  LOG_DEBUG<<"check simu table :" << endm;
305  for (i = 0 ; i < mSsdSimuData->GetNRows(); i++)
306  {
307  nStrip = (int)(strip[i].id_strip/100000.);
308  idWaf = strip[i].id_strip-10000*((int)(strip[i].id_strip/10000.));
309  iWaf = idWaferToWafer(idWaf);
310  iLad = (int)(idWaf - mSsdLayer*1000 - (iWaf+1)*100 - 1);
311  iSide = (strip[i].id_strip - nStrip*100000 - idWaf)/10000;
312  LOG_DEBUG<<Form("side=%d idWafer=%d Ladder=%d wafer=%d nstrip=%d signal=%d",iSide,idWaf,iLad,iWaf,nStrip,strip[i].adc_count)<<endm;
313  for(Int_t ii=0;ii<5;ii++){printf("id_mchit[%d] =%d\n",ii,strip[i].id_mchit[ii]);}
314  }
315  if(Debug()){
316  LOG_DEBUG << "check real table : " << endm;
317  strip = mSsdrealData->GetTable();
318  for (i = 0 ; i < mSsdrealData->GetNRows(); i++)
319  {
320  nStrip = (int)(strip[i].id_strip/100000.);
321  idWaf = strip[i].id_strip-10000*((int)(strip[i].id_strip/10000.));
322  iWaf = idWaferToWafer(idWaf);
323  iLad = (int)(idWaf - mSsdLayer*1000 - (iWaf+1)*100 - 1);
324  iSide = (strip[i].id_strip - nStrip*100000 - idWaf)/10000;
325  LOG_DEBUG<<Form("side=%d Ladder=%d wafer=%d nstrip=%d signal=%d",iSide,iLad,iWaf,nStrip,strip[i].adc_count)<<endm;
326  }
327  }
328 }
329 //__________________________________________________________________________
330 void StSsdEmbeddingMaker::setSsdParameters(ssdDimensions_st *geom_par){
331  mDimensions = geom_par;
332  mSsdLayer = 7;
333  mDetectorLargeEdge = 2.*geom_par[0].waferHalfActLength;
334  mDetectorSmallEdge = 2.*geom_par[0].waferHalfActWidth;
335  mNLadder = 20;
336  mNWaferPerLadder = geom_par[0].wafersPerLadder;
337  mNStripPerSide = geom_par[0].stripPerSide;
338  mStripPitch = geom_par[0].stripPitch;
339  mTheta = geom_par[0].stereoAngle;
340 }
341 //___________________________________________________________________________
342 
TDataSet::GetObject
virtual TObject * GetObject() const
The depricated method (left here for the sake of the backward compatibility)
Definition: TDataSet.cxx:428
kStOK
Definition: Stypes.h:40
kStErr
Definition: Stypes.h:44
kStOk
Definition: Stypes.h:41
StSsdEmbeddingMaker::InitRun
virtual Int_t InitRun(int)
All database dependent data are read here.
Definition: StSsdEmbeddingMaker.cxx:60
Generated by   doxygen 1.8.5