StRoot  1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
StFcsEventDisplay.cxx
1 // \class StFmsEventDisplay
2 // \author Akio Ogawa
3 //
4 // $Id: StFcsEventDisplay.cxx,v 1.1 2021/03/30 13:34:15 akio Exp $
5 // $Log: StFcsEventDisplay.cxx,v $
6 // Revision 1.1 2021/03/30 13:34:15 akio
7 // Moved from $CVSROOT/offline/upgrade/akio/ to $CVSROOT/StRoot/StSpinPool/
8 //
9 // Revision 1.15 2021/02/25 21:27:50 akio
10 // Int_t -> int
11 //
12 // Revision 1.14 2021/02/13 21:37:31 akio
13 // #ifdef ___USESTGC___
14 //
15 // Revision 1.13 2020/05/29 18:54:29 akio
16 // Adding EPD as PRES
17 //
18 // Revision 1.12 2019/10/23 19:20:50 akio
19 // forgot to add display offset for point
20 //
21 // Revision 1.11 2019/10/23 13:30:49 akio
22 // adding display of FcsPoint
23 //
24 // Revision 1.10 2019/08/01 18:38:22 akio
25 // Added STGC
26 //
27 // Revision 1.9 2019/07/02 14:43:16 akio
28 // fixed TColor id ovewrite problem, fixed pdf file creation
29 //
30 // Revision 1.8 2019/06/27 16:13:27 akio
31 // minor updates on not printing energy for ecal (unless run19) for visibility
32 //
33 // Revision 1.7 2019/06/26 18:01:38 akio
34 // make clusterid/energy printing smaller for run21
35 //
36 // Revision 1.6 2019/06/25 16:39:30 akio
37 // Added run19 version, and drawing cluster Id and tower energy
38 //
39 // Revision 1.5 2019/06/21 18:53:25 akio
40 // added run19 version
41 //
42 // Revision 1.4 2019/06/07 18:13:23 akio
43 // added filter
44 //
45 // Revision 1.3 2019/05/17 15:58:06 akio
46 // removing TApplication
47 //
48 // Revision 1.2 2019/05/16 17:24:01 akio
49 // minor bug fix drawing det=6 which doesn't exist
50 //
51 // Revision 1.1 2018/11/14 16:50:14 akio
52 // FCS codes in offline/upgrade/akio
53 //
54 
55 #include "StFcsEventDisplay.h"
56 
57 #include "StEvent/StEnumerations.h"
58 #include "StMessMgr.h"
59 #include "Stypes.h"
60 #include "StEventTypes.h"
61 #include "StMuDSTMaker/COMMON/StMuTypes.hh"
62 
63 #include "StThreeVectorF.hh"
64 #include "StFcsDbMaker/StFcsDb.h"
65 #include "StRoot/StEpdUtil/StEpdGeom.h"
66 
67 #ifdef ___USESTGC___
68 #include "StStgcDbMaker/StStgcDbMaker.h"
69 #endif
70 
71 //#include "TApplication.h"
72 #include "TCanvas.h"
73 #include "TFile.h"
74 #include "TH1F.h"
75 #include "TH2F.h"
76 #include "TLine.h"
77 #include "TBox.h"
78 #include "TMarker.h"
79 #include "TString.h"
80 #include "TColor.h"
81 #include "TStyle.h"
82 #include "TText.h"
83 #include "TROOT.h"
84 
85 ClassImp(StFcsEventDisplay)
86 
87 StFcsEventDisplay::StFcsEventDisplay(const char* name):
88  StMaker(name),mFilename((char *)"fcsEventDisplay.pdf"){}
89 
90 StFcsEventDisplay::~StFcsEventDisplay(){}
91 
92 int StFcsEventDisplay::Init(){
93  mFcsDb=static_cast<StFcsDb*>(GetDataSet("fcsDb"));
94  if(!mFcsDb){
95  LOG_ERROR << "StFcsEventDisplay::InitRun Failed to get StFcsDb" << endm;
96  return kStFatal;
97  }
98 
99 #ifdef ___USESTGC___
100  mStgcDbMaker=static_cast<StStgcDbMaker*>(GetMaker("stgcDb"));
101  if(!mStgcDbMaker){
102  LOG_WARN << "StFcsEventDisplay::InitRun Failed to get StStgcDbMaker" << endm;
103  //return kStFatal;
104  }
105 #endif
106 
107  //mApplication = new TApplication("EvtDsp",0,0);
108  mCanvas=new TCanvas("FCSEventtDisplay","FCSEventtDisplay",10,10,2000,2000);
109  gStyle->SetOptStat(0);
110 
111  //EPD geom
112  mEpdgeo=new StEpdGeom;
113 
114  return kStOK;
115 }
116 
118  //mApplication->Terminate();
119  return kStOK;
120 }
121 
122 void setColor(int id, float v, float alpha=1.0){ //blue->green for 0->0.5 and green->red for 0.5->1.0
123  static const float PI=3.141592654;
124  TColor *color = gROOT->GetColor(id);
125  if(!color) color=new TColor(id,0,0,0);
126  float r=0.0,g=0.0,b=0.0;
127  if(v>1.0) v=1.0;
128  if(v<0.0) v=0.0;
129  if(v>0.5){
130  r=sin((v-0.5)/0.5*PI/2.0);
131  g=cos((v-0.5)/0.5*PI/2.0);
132  }else{
133  g=sin(v/0.5*PI/2.0);
134  b=cos(v/0.5*PI/2.0);
135  }
136  color->SetRGB(r,g,b);
137  color->SetAlpha(alpha);
138  //printf("%d %f %f %f %f\n",id,v,r,g,b);
139 }
140 
141 void scale(float x=160.0, float dx=5.0, float ymin=-210.0, float ymax=210.0, int ny=40){
142  for(int i=0; i<ny; i++){
143  setColor(5000+i,float(i+0.5)/float(ny));
144  float dy=(ymax-ymin)/float(ny);
145  TBox* b=new TBox(x,ymin+dy*i,x+dx,ymin+dy*(i+1));
146  b->SetFillColor(5000+i);
147  b->Draw();
148  }
149 }
150 
152  StEvent* event = (StEvent*)GetInputDS("StEvent");
153  if(!event) {LOG_ERROR << "StFcsEventDisplay::Make did not find StEvent"<<endm; return kStErr;}
154  mFcsColl = event->fcsCollection();
155  if(!mFcsColl) {LOG_ERROR << "StFcsEventDisplay::Make did not find StEvent->StFcsCollection"<<endm; return kStErr;}
156 
157 #ifdef ___USESTGC___
158  mStgcColl = event->stgcCollection();
159  if(!mStgcColl) {LOG_ERROR << "StFcsEventDisplay::Make did not find StEvent->StStgcCollection"<<endm;}
160 #endif
161 
162  if(mNAccepted < mMaxEvents){
163  mNEvents++;
164  if(mFilter==1 &&
165  mFcsColl->numberOfHits(0)+ mFcsColl->numberOfHits(1)+ mFcsColl->numberOfHits(2)+ mFcsColl->numberOfHits(3)==0) return kStOK;
166  mNAccepted++;
167 
168  mCanvas->Clear();
169  char cc[100];
170  sprintf(cc,"FCSEventDisplayEvt=%d",mNEvents);
171  mFcsColl->print();
172 
173  //frame
174  float xmin=-160;
175  float xmax=160;
176  float ymin=-220;
177  float ymax=220;
178  float yoff[kFcsNDet]={100.0,100.0,-100.0,-100.0,100.0,100.0};
179  float xoff[kFcsNDet]={ 0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
180  float stgcOffX=0.0;
181  float stgcOffY=0.0;
182  float cTxSize=0.01;
183  float eTxSize=0.007;
184  TH2F* frame;
185  if(mRun19>0){
186  xmin=0;
187  xmax=120;
188  ymin=-110;
189  ymax=10;
190  xoff[1]=0;
191  xoff[3]=-50;
192  xoff[5]=0;
193  yoff[1]=0;
194  yoff[3]=0;
195  yoff[5]=40;
196  cTxSize=0.02;
197  eTxSize=0.015;
198  }
199  frame=new TH2F(cc,cc,1,xmin,xmax,1,ymin,ymax);
200  float dy=ymax-ymin;
201 
202  frame->Draw();
203  scale(xmax,5,ymin,ymax);
204  if(mMaxE==2 && mMinE==-1){
205  TText* t1=new TText(xmax+5, ymax, "100GeV"); t1->SetTextSize(0.02); t1->Draw();
206  TText* t2=new TText(xmax+5, ymin+dy*2/3.0, "10GeV" ); t2->SetTextSize(0.02); t2->Draw();
207  TText* t3=new TText(xmax+5, ymin+dy/3.0, "1GeV" ); t3->SetTextSize(0.02); t3->Draw();
208  TText* t4=new TText(xmax+5, ymin, "0.1GeV"); t4->SetTextSize(0.02); t4->Draw();
209  }else if(mMaxE==1 && mMinE==-3){
210  TText* t1=new TText(xmax+5, ymax, "1GeV"); t1->SetTextSize(0.02); t1->Draw();
211  TText* t2=new TText(xmax+5, ymin+dy*2/3.0, "100MeV" ); t2->SetTextSize(0.02); t2->Draw();
212  TText* t3=new TText(xmax+5, ymin+dy/3.0, "10MeV" ); t3->SetTextSize(0.02); t3->Draw();
213  TText* t4=new TText(xmax+5, ymin, "1MeV"); t4->SetTextSize(0.02); t4->Draw();
214  }else{
215  TText* t1=new TText(xmax+5, ymax, Form("1E%2.0f GeV",mMaxE)); t1->SetTextSize(0.02); t1->Draw();
216  TText* t4=new TText(xmax+5, ymin, Form("1E%2.0f GeV",mMinE)); t4->SetTextSize(0.02); t4->Draw();
217  }
218 
219  int color = 1000;
220  for(int det=0; det<kFcsNDet; det++){
221  StSPtrVecFcsHit& hits = mFcsColl->hits(det);
222  StSPtrVecFcsCluster& clusters = mFcsColl->clusters(det);
223  StSPtrVecFcsPoint& points = mFcsColl->points(det);
224  int nh=mFcsColl->numberOfHits(det);
225  int nc=mFcsColl->numberOfClusters(det);
226  int np=mFcsColl->numberOfPoints(det);
227  if(mDebug>0) LOG_INFO << Form("StFcsEventDisplay Det=%1d nhit=%4d nclu=%3d",det,nh,nc) << endm;
228 
229  //first all towers
230  if(det<=kFcsHcalSouthDetId){ //ecal & hcal
231  for(int id=0; id<mFcsDb->maxId(det); id++){
232  StThreeVectorF xyz = mFcsDb->getStarXYZ(det,id);
233  float wx=mFcsDb->getXWidth(det);
234  float wy=mFcsDb->getYWidth(det);
235  TBox* cell=new TBox(xyz.x()-wx/2.0+xoff[det], xyz.y()+yoff[det]-wy/2.0,
236  xyz.x()+wx/2.0+xoff[det], xyz.y()+yoff[det]+wy/2.0);
237  cell->SetFillStyle(0);
238  //cell->SetFillColorAlpha(0,0.9);
239  cell->SetLineColor(1); cell->SetLineWidth(1);
240  cell->Draw("l");
241  }
242  }
243 
244  //draw hits
245  for(int i=0; i<nh; i++){
246  StFcsHit* hit=hits[i];
247  if(det<=kFcsHcalSouthDetId){ //ecal & hcal
248  StThreeVectorF xyz = mFcsDb->getStarXYZ(hit);
249  float wx=mFcsDb->getXWidth(det);
250  float wy=mFcsDb->getYWidth(det);
251  TBox* cell=new TBox(xyz.x()-wx/2.0+xoff[det], xyz.y()+yoff[det]-wy/2.0,
252  xyz.x()+wx/2.0+xoff[det], xyz.y()+yoff[det]+wy/2.0);
253  float e=hit->energy();
254  if(e<=0.0) continue;
255  float logE=log10(e);
256  if(logE>mMaxE) logE=mMaxE;
257  if(logE<mMinE) logE=mMinE;
258  setColor(color,(logE - mMinE)/(mMaxE-mMinE));
259  cell->SetFillColor(color);
260  cell->Draw();
261  color++;
262  if(mRun19 || det>=kFcsHcalNorthDetId){
263  TText* tx=new TText(xyz.x()+xoff[det],xyz.y()+yoff[det]-wy/4.0,Form("%4.2f",e));
264  tx->SetTextSize(eTxSize); tx->SetTextAlign(22);
265  tx->Draw();
266  }
267  }else if(det==kFcsPresNorthDetId || det==kFcsPresSouthDetId){//EPD as Pres
268  double zepd=375.0;
269  double zfcs=710.0+13.90+15.0;
270  double zr=zfcs/zepd;
271  int pp,tt,id,n;
272  double x[5],y[5];
273  mFcsDb->getEPDfromId(det,hit->id(),pp,tt);
274  mEpdgeo->GetCorners(100*pp+tt,&n,x,y);
275  for(int i=0; i<n; i++){
276  int j=i+1;
277  if(i==n-1) j=0;
278  TLine *l=new TLine(zr*x[i]+xoff[det],zr*y[i]+yoff[det],
279  zr*x[j]+xoff[det],zr*y[j]+yoff[det]);
280  l->SetLineColor(kRed);
281  l->SetLineWidth(4);
282  l->Draw();
283  }
284  }
285  if(mDebug>0) LOG_INFO << Form("StFcsEventDisplay Det=%1d nhit=%4d i=%3d color=%4d",det,nh,i,color) << endm;
286  }
287 
288  //clusters
289  for(int i=0; i<nc; i++){
290  StFcsCluster* clu=clusters[i];
291  float nt=float(clu->nTowers());
292  for(int j=0; j<nt; j++){
293  StFcsHit* hit = clu->hits()[j];
294  StThreeVectorF xyz = mFcsDb->getStarXYZ(hit);
295  float wx=mFcsDb->getXWidth(det);
296  float wy=mFcsDb->getYWidth(det);
297  TBox* cell=new TBox(xyz.x()-wx/2.0+xoff[det], xyz.y()+yoff[det]-wy/2.0,
298  xyz.x()+wx/2.0+xoff[det], xyz.y()+yoff[det]+wy/2.0);
299  cell->SetFillStyle(0);
300  cell->SetLineColor(2+i);
301  cell->SetLineWidth(2);
302  cell->Draw("l");
303  TText* tx;
304  if(mRun19 || det>=kFcsHcalNorthDetId){
305  tx=new TText(xyz.x()+xoff[det],xyz.y()+yoff[det]+wy/4.0,Form("%d",i));
306  }else{
307  tx=new TText(xyz.x()+xoff[det],xyz.y()+yoff[det],Form("%d",i));
308  }
309  tx->SetTextSize(cTxSize); tx->SetTextAlign(22);
310  tx->Draw();
311  }
312  StThreeVectorF cxyz=mFcsDb->getStarXYZfromColumnRow(det,clu->x(),clu->y());
313  TMarker* m=new TMarker(cxyz.x()+xoff[det],cxyz.y()+yoff[det],26);
314  m->SetMarkerColor(kMagenta);
315  m->SetMarkerSize(2);
316  m->Draw();
317  TText* tx;
318  tx=new TText(cxyz.x()+xoff[det],cxyz.y()+yoff[det],
319  Form(" %d/%3.1f/%3.1f/%3.0f",
320  clu->nTowers(),clu->sigmaMax(),clu->sigmaMin(),
321  clu->theta()*180/3.141592654));
322  tx->SetTextSize(cTxSize);
323  tx->SetTextColor(kBlue);
324  tx->Draw();
325  int nphoton=clu->nPoints();
326  for(int j=0; j<nphoton; j++){
327  StFcsPoint* point = clu->points()[j];
328  StThreeVectorF xyz=point->xyz();
329  TMarker* m=new TMarker(xyz.x()+xoff[det],xyz.y()+yoff[det],29);
330  m->SetMarkerSize(2);
331  m->SetMarkerColor(kRed);
332  m->Draw();
333  }
334  }
335  }
336 
337 #ifdef ___USESTGC___
338  //STGC
339  if(mStgcColl){
340  for(int det=0; det<kStgcNDet; det++){
341  StSPtrVecFtsStgcHit& hits = mStgcColl->hits(det);
342  int ns=mStgcColl->numberOfHits(det);
343  if(mDebug>0) LOG_INFO << Form("StFcsEventDisplay STGC det=%1d nhit=%4d",det,ns) << endm;
344  for(int i=0; i<ns; i++){
345  StFtsStgcHit* hit=hits[i];
346  unsigned short fee = hit->fee();
347  unsigned short altro= hit->altro();
348  unsigned short ch = hit->ch();
349  unsigned int id = mStgcDbMaker->toId(fee,altro,ch);
350  float x,y,dx,dy,z,dz;
351  mStgcDbMaker->globalPosition(id,x,y,dx,dy,z,dz);
352  TBox* cell=new TBox(x-dx/2+stgcOffX, y-dy/2+stgcOffY,
353  x+dx/2+stgcOffX, y+dy/2+stgcOffY);
354  cell->SetFillStyle(0);
355  cell->SetLineColor(kRed);
356  cell->SetLineWidth(1);
357  cell->Draw();
358  }
359  }
360 
361  }
362 #endif
363 
364  mCanvas->Update();
365  TString f(mFilename);
366  if(f.Contains(".pdf")){ //not working
367  if(mNAccepted==1) f.Append("(");
368  else if(mNAccepted==mMaxEvents) f.Append(")");
369  LOG_INFO << "Saving " << f.Data() << endm;
370  mCanvas->Print(f.Data(),"pdf");
371  }else if(f.Contains(".png")){
372  f.ReplaceAll("eventDisplay.png",Form("%d.eventDisplay.png",mNEvents));
373  LOG_INFO << "Saving " << f.Data() << endm;
374  mCanvas->SaveAs(f.Data());
375  }
376  }
377  return kStOK;
378 }
StThreeVectorD getStarXYZfromColumnRow(int det, float col, float row, float FcsZ=-1.0) const
get the STAR frame cooridnates from other way
Definition: StFcsDb.cxx:863
float getXWidth(int det) const
get the angle of the detector
Definition: StFcsDb.cxx:627
StThreeVectorD getStarXYZ(int det, float FcsX, float FcsY, float FcsZ=-1.0, float zVertex=0.0) const
get the STAR frame cooridnates from local XYZ [cm]
Definition: StFcsDb.cxx:868
void getEPDfromId(int det, int id, int &pp, int &tt)
Get FCS&#39;s EPD map foom EPD mapping.
Definition: StFcsDb.cxx:1484
Definition: Stypes.h:40
int maxId(int det) const
number of column
Definition: StFcsDb.cxx:468
Definition: Stypes.h:44
float getYWidth(int det) const
get the X width of the cell
Definition: StFcsDb.cxx:635