StRoot  1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
StFcsEcalHcalMipMaker.cxx
1 #include "StFcsEcalHcalMipMaker.h"
2 #include "TDataSetIter.h"
3 #include "StDAQMaker/StDAQReader.h"
4 
5 #include "StRoot/StEvent/StEvent.h"
6 #include "StRoot/St_base/StMessMgr.h"
7 #include "StRoot/StEvent/StTriggerData.h"
8 #include "StRoot/StEvent/StFcsCollection.h"
9 #include "StRoot/StEvent/StFcsHit.h"
10 #include "StRoot/StEvent/StFcsCluster.h"
11 #include "StRoot/StFcsDbMaker/StFcsDb.h"
12 #include "StRoot/StSpinPool/StFcsRawDaqReader/StFcsRawDaqReader.h"
13 
14 #include "TH1F.h"
15 #include "TH2F.h"
16 #include "TString.h"
17 #include "TFile.h"
18 #include "TCanvas.h"
19 
20 #include <string.h>
21 #include <time.h>
22 
23 //ClassImp(StFcsEcalHcalMipMaker)
24 
25 //_____________________________________________________________________________
27 
29 
30 }
31 
32 //_____________________________________________________________________________
34 
36 
37 
38 //_____________________________________________________________________________
41  mFcsDb = static_cast<StFcsDb*>(GetDataSet("fcsDb"));
42 
43  sprintf(mFilename,"%8d.mip.root",mRun);
44  LOG_INFO << "StFcsEcalHcalMipMaker::Init - Opening "<<mFilename<<endm;
45  mFile=new TFile(mFilename,"RECREATE");
46 
47  const char* nameEH[2] = {"Ecal","Hcal"};
48  const char* nameNS[kFcsNorthSouth] = {"N","S"};
49  char f[100], t[100]; //histogram file names
50 
51  for(int eh=0; eh<2; eh++){
52  mNClu[eh] = new TH1F(Form("NCluster_%s",nameEH[eh]),
53  Form("NCluster_%s",nameEH[eh]),
54  50,0,50);
55  mNTowClu[eh] = new TH1F(Form("NTowerCluster_%s",nameEH[eh]),
56  Form("NTowerCluster_%s",nameEH[eh]),
57  30,0,30);
58  mNNeiClu[eh] = new TH1F(Form("NNeiCluster_%s",nameEH[eh]),
59  Form("NNeiCluster_%s",nameEH[eh]),
60  20,0,20);
61  }
62  const int nbin=256;
63  const float max=3.0;
64  for(int det=0; det<4; det++){
65  int ns = mFcsDb->northSouth(det);
66  int eh = mFcsDb->ecalHcalPres(det);
67  int maxid = mFcsDb->maxId(det);
68  mAdc[det] = new TH2F(Form("ADC_%s%s",nameEH[eh],nameNS[ns]),
69  Form("ADC %s%s",nameEH[eh],nameNS[ns]),
70  maxid,0,maxid,nbin,0,max);
71  mAdcSingleTower[det] = new TH2F(Form("ADCSingle_%s%s",nameEH[eh],nameNS[ns]),
72  Form("ADCSingle %s%s",nameEH[eh],nameNS[ns]),
73  maxid,0,maxid,nbin,0,max);
74  mAdcEcalMatch[det] = new TH2F(Form("ADCEcal_%s%s",nameEH[eh],nameNS[ns]),
75  Form("ADCEcalCoin %s%s",nameEH[eh],nameNS[ns]),
76  maxid,0,maxid,nbin,0,max);
77  }
78  mX = new TH2F("X","X; XEcal; XHcal",100,-100,100,100,-100,100);
79  mY = new TH2F("Y","Y; YEcal; YHcal",100, 0,200,100, 0,200);
80  mDXX= new TH2F("DXX","DX; XHcal; DX",100,-100,100,100,-30,30);
81  mDXY= new TH2F("DXY","DX; YHcal; DX",100, 0,200,100,-30,30);
82  mDYX= new TH2F("DYX","DY; XHcal; DY",100,-100,100,100,-30,30);
83  mDYY= new TH2F("DYY","DY; YHcal; DY",100, 0,200,100,-30,30);
84  mDR= new TH2F("DR","DR; RHcal; DR",100,0,100,100,0,30);
85  return kStOk;
86 }
87 
88 //_____________________________________________________________________________
91  mFcsCollection=0;
92  StEvent* event= (StEvent*)GetInputDS("StEvent");
93  if(!event) {
94  LOG_INFO << "No StEvent found" << endm;
95  }else{
96  mFcsCollection=event->fcsCollection();
97  }
98  if(!mFcsCollection){
99  LOG_INFO << "No StFcsCollection found" << endm;
100  return kStErr;
101  }
102 
103  const int maxTower=1;
104  const float ecalELow=0.15;
105  const float ecalEHigh=0.50;
106  const float dRCut=15.0;
107 
108  for(int ns=0; ns<2; ns++){
109  StSPtrVecFcsCluster& ecal= mFcsCollection->clusters(ns);
110  StSPtrVecFcsCluster& hcal= mFcsCollection->clusters(ns+2);
111  int nEcal=mFcsCollection->numberOfClusters(ns);
112  int nHcal=mFcsCollection->numberOfClusters(ns+2);
113  mNClu[0]->Fill(nEcal);
114  mNClu[1]->Fill(nHcal);
115 
116  //loop over Ecal clusters
117  for(int i=0; i<nEcal;i++){
118  StFcsCluster* clu=ecal[i];
119  int id=clu->hits()[0]->id();
120  int ntow=clu->nTowers();
121  int nnei=clu->nNeighbor();
122  mNTowClu[0]->Fill(ntow);
123  mNNeiClu[0]->Fill(nnei);
124  mAdc[ns]->Fill(id,clu->energy());
125  if(ntow<=maxTower && nnei==0){
126  mAdcSingleTower[ns]->Fill(id,clu->energy());
127  }
128  }
129 
130  //loop over Hcal clusters
131  for(int i=0; i<nHcal;i++){
132  StFcsCluster* clu=hcal[i];
133  int id=clu->hits()[0]->id();
134  int ntow=clu->nTowers();
135  int nnei=clu->nNeighbor();
136  mNTowClu[1]->Fill(ntow);
137  mNNeiClu[1]->Fill(nnei);
138  mAdc[ns+2]->Fill(id,clu->energy());
139  if(ntow<=2 && nnei==0){
140  mAdcSingleTower[ns+2]->Fill(id,clu->energy());
141  for(int i=0; i<nEcal;i++){ //loop over Ecal
142  StFcsCluster* eclu=ecal[i];
143  int entow=eclu->nTowers();
144  int ennei=eclu->nNeighbor();
145  if(entow<=maxTower && ennei==0 && eclu->energy()>ecalELow && eclu->energy()<ecalEHigh){
146  double px = mFcsDb->getHcalProjectedToEcalX(ns,clu->x()) * (ns*2-1);
147  double py = mFcsDb->getHcalProjectedToEcalY(ns,clu->y());
148  double ex = eclu->x() * mFcsDb->getXWidth(ns) * (ns*2-1);
149  double ey = eclu->y() * mFcsDb->getYWidth(ns);
150  double dx=ex-px;
151  double dy=ey-py;
152  double pr=sqrt(px*px+py*py);
153  double dr=sqrt(dx*dx+dy*dy);
154  mX->Fill(ex,px);
155  mY->Fill(ey,py);
156  if(abs(dy)<dRCut){
157  mDXX->Fill(px,dx);
158  mDXY->Fill(py,dx);
159  }
160  if(abs(dx)<dRCut){
161  mDYX->Fill(px,dy);
162  mDYY->Fill(py,dy);
163  }
164  mDR->Fill(pr,dr);
165  if(dr<dRCut){
166  mAdcEcalMatch[ns+2]->Fill(id,clu->energy());
167  }
168  }
169  }
170  }
171  }
172  }
173  return kStOK;
174 }
175 
177  mFile->Write();
178  mFile->Close();
179  printf("StFcsEcalHcalMipMaker::Finish - Closing %s\n",mFilename);
180  return kStOK;
181 };
182 
183 ClassImp(StFcsEcalHcalMipMaker);
184 
virtual Int_t Init()
Init - the top level StChain calls to initialize all its makers.
virtual Int_t Make()
Make - this method is called in loop for each event.
int northSouth(int det) const
Ecal=0, Hcal=1, Pres=2.
Definition: StFcsDb.cxx:430
double getHcalProjectedToEcalX(int ns, double hcalLocalX, double zvtx=0.0)
Project Hcal local X to Ecal local X [cm].
Definition: StFcsDb.cxx:695
StFcsEcalHcalMipMaker(const char *name="FcsEcalHcalMip")
constructor
float getXWidth(int det) const
get the angle of the detector
Definition: StFcsDb.cxx:627
int ecalHcalPres(int det) const
Ecal North=0, Ecal South=1, Hcal North=2, Hcal South=3, Pres=4/5.
Definition: StFcsDb.cxx:423
Definition: Stypes.h:40
int maxId(int det) const
number of column
Definition: StFcsDb.cxx:468
Definition: Stypes.h:44
virtual ~StFcsEcalHcalMipMaker()
destructor
Definition: Stypes.h:41
float getYWidth(int det) const
get the X width of the cell
Definition: StFcsDb.cxx:635
double getHcalProjectedToEcalY(int ns, double hcalLocalY, double zvtx=0.0)
Project Hcal local Y to Ecal local Y [cm].
Definition: StFcsDb.cxx:713