er  dev
ERHe8EventHeader.cxx
1 #include <iostream>
2 
3 #include "ERHe8EventHeader.h"
4 
5 ERHe8EventHeader::ERHe8EventHeader(){
6 
7 }
8 //-----------------------------------------------------------------------------
9 Bool_t ERHe8EventHeader::Register(TTree* tree, TString branchName){
10  tree->SetMakeClass(1);
11  tree->SetBranchAddress(branchName + TString(".nevent"),&this->HE8Event_nevent);
12  tree->SetBranchAddress(branchName + TString(".trigger"),&HE8Event_trigger);
13  tree->SetBranchAddress(branchName + TString(".subevents"),&HE8Event_subevents);
14  tree->SetBranchAddress(branchName + TString(".evsize"),&HE8Event_evsize);
15  tree->SetBranchAddress(branchName + TString(".time"),&HE8Event_time);
16 
17  ReadInputFile();
18  ReactionPreparation();
19 
20  return kTRUE;
21 }
22 //-----------------------------------------------------------------------------
23 void ERHe8EventHeader::ReadInputFile(){
24  std::cout << "Filling header from input files" << std::endl;
25  char Zeros[32];
26  ReIN.Simulation = false;
27  ReIN.Vertex = false;
28  ReIN.DetectorTune = false;
29  ReIN.ToFis = false;
30  ReIN.TRACKINGis = false;
31  ReIN.TrackCheck = false;
32  ReIN.WriteRawData = false;
33  ReIN.WriteCalData = false;
34  ReIN.WriteToFData = false;
35  ReIN.WriteTrackData = false;
36  ReIN.WriteTelData = false;
37  ReIN.WritePhysData = false;
38  ReIN.WritePlayData = false;
39  ReIN.WriteRawTrack = false;
40  /******************** Readout ReactionInput.dat ************************/
41  printf("************************************************************\n");
42  TString filePath = gSystem->Getenv("VMCWORKDIR") + TString("/input/") + fReactionInputFile;
43  FILE *F1 = fopen(filePath.Data(),"r");
44  if(F1==NULL) {
45  printf("ERHe8EventHeader: File ReactionInput.dat was not found\n");
46  } else {
47  int res;
48  res = fscanf(F1,"%s\n",ReIN.ReactionName);
49  res = fscanf(F1,"%s %s\n",ReIN.AboutBeam,ReIN.EnergyUn);
50  res = fscanf(F1,"%s %s\n",ReIN.AboutSlit,ReIN.SlitUn);
51  res = fscanf(F1,"%s %s\n",Zeros,ReIN.Mechanism);
52  res = fscanf(F1,"%s %s\n",Zeros,ReIN.Fname);
53  res = fscanf(F1,"%s %s\n",Zeros,Zeros);
54  if(!strcmp(Zeros,"yes")) {ReIN.Simulation = true;}
55  res = fscanf(F1,"%s %s\n",Zeros,Zeros);
56  if(!strcmp(Zeros,"yes")) {ReIN.Vertex = true;}
57  res = fscanf(F1,"%s %s\n",Zeros,Zeros);
58  if(!strcmp(Zeros,"yes")) {ReIN.DetectorTune = true;}
59  res = fscanf(F1,"%s %s\n",Zeros,Zeros);
60  if(!strcmp(Zeros,"yes")) {ReIN.ToFis = true;}
61  res = fscanf(F1,"%s %s\n",Zeros,Zeros);
62  if(!strcmp(Zeros,"yes")) {ReIN.TRACKINGis = true;}
63  res = fscanf(F1,"%s %s\n",Zeros,Zeros);
64  if(!strcmp(Zeros,"yes")) {ReIN.TrackCheck = true;}
65  res = fscanf(F1,"%s %s\n",Zeros,Zeros);
66  if(!strcmp(Zeros,"yes")) {ReIN.WriteRawData = true;}
67  res = fscanf(F1,"%s %s\n",Zeros,Zeros);
68  if(!strcmp(Zeros,"yes")) {ReIN.WriteCalData = true;}
69  res = fscanf(F1,"%s %s\n",Zeros,Zeros);
70  if(!strcmp(Zeros,"yes")) {ReIN.WriteToFData = true;}
71  res = fscanf(F1,"%s %s\n",Zeros,Zeros);
72  if(!strcmp(Zeros,"yes")) {ReIN.WriteTrackData = true;}
73  res = fscanf(F1,"%s %s\n",Zeros,Zeros);
74  if(!strcmp(Zeros,"yes")) {ReIN.WriteTelData = true;}
75  res = fscanf(F1,"%s %s\n",Zeros,Zeros);
76  if(!strcmp(Zeros,"yes")) {ReIN.WritePhysData = true;}
77  res = fscanf(F1,"%s %s\n",Zeros,Zeros);
78  if(!strcmp(Zeros,"yes")) {ReIN.WritePlayData = true;}
79  res = fscanf(F1,"%s %s\n",Zeros,Zeros);
80  if(!strcmp(Zeros,"yes")) {ReIN.WriteRawTrack = true;}
81 
82  for(int m=0;m<8;m++)
83  {
84  int res2;
85  res2 = fscanf(F1,"%s %i\n",Zeros,&ReIN.ifill[m]);
86  }
87  printf("Main: File ReactionInput.dat has been read\n");
88  }
89  fclose(F1);
90 
91  /*********************** Readout MWPC parameters:***********************/
92  if(ReIN.TRACKINGis)
93  {
94  filePath = gSystem->Getenv("VMCWORKDIR") + TString("/input/track.dat");
95  FILE *F2 = fopen(filePath.Data(),"r");
96  if(F2==NULL)
97  {
98  printf("Main: File track.dat was not found\n");
99  } else {
100  int res;
101  res = fscanf(F1,"%s %s %lf\n",Zeros,UpMat.MWwinMatter,&UpMat.MWwinThick);
102  res = fscanf(F1,"%s %s %lf\n",Zeros,UpMat.MWgasMatter,&UpMat.MWgasThick);
103  res = fscanf(F1,"%s %s %lf\n",Zeros,UpMat.MWcathMatter,&UpMat.MWcathThick);
104  res = fscanf(F1,"%s %i\n",Zeros,&UpMat.MWNcathodes);
105  res = fscanf(F1,"%s %lf\n",Zeros,&UpMat.MWXYdist);
106  res = fscanf(F1,"%s %lf\n",Zeros,&UpMat.MWstep);
107  res = fscanf(F1,"%s %i\n",Zeros,&UpMat.MWNwires);
108  res = fscanf(F1,"%s\n",Zeros);
109  res = fscanf(F1,"%s %lf\n",Zeros,&UpMat.MWfarDist);
110  res = fscanf(F1,"%s %lf %lf\n",Zeros,&UpMat.MWfarXshift,&UpMat.MWfarYshift);
111  res = fscanf(F1,"%s %i %i\n",Zeros,&UpMat.MWfarXNum,&UpMat.MWfarYNum);
112  res = fscanf(F1,"%s\n",Zeros);
113  res = fscanf(F1,"%s %lf\n",Zeros,&UpMat.MWclosDist);
114  res = fscanf(F1,"%s %lf %lf\n",Zeros,&UpMat.MWclosXshift,&UpMat.MWclosYshift);
115  res = fscanf(F1,"%s %i %i\n",Zeros,&UpMat.MWclosXNum,&UpMat.MWclosYNum);
116  printf("Main: File track.dat has been read\n");
117  }
118  fclose(F2);
119  }
120  UpMat.MWcathThick *= UpMat.MWNcathodes;
121 
122  //TODO Check that this section works correctly!!! ========================
123  /*********************** Readout ToF parameters:************************/
124  if(ReIN.ToFis)
125  {
126  filePath = gSystem->Getenv("VMCWORKDIR") + TString("/input/tof.dat");
127  FILE *F3 = fopen(filePath.Data(),"r");
128  if(F3==NULL) {
129  printf("Main: File tof.dat was not found\n");
130  }
131  else
132  {
133  int res;
134  res = fscanf(F1,"%s %lf %lf\n",UpMat.PlasticMatter1,&UpMat.PlasticAngle1,&UpMat.PlasticThick1);
135  res = fscanf(F1,"%s %lf %lf\n",UpMat.PlasticMatter2,&UpMat.PlasticAngle2,&UpMat.PlasticThick2);
136  res = fscanf(F1,"%s %lf\n",Zeros,&UpMat.PlasticDist);
137  res = fscanf(F1,"%s %lf\n",Zeros,&UpMat.ToFRes);
138  res = fscanf(F1,"%s %lf %lf\n",Zeros,&UpMat.tF3l_rng,&UpMat.tF3r_rng);
139  res = fscanf(F1,"%s %lf %lf\n",Zeros,&UpMat.tF4l_rng,&UpMat.tF4r_rng);
140  res = fscanf(F1,"%s %lf %lf\n",Zeros,&UpMat.tF3_dlt,&UpMat.tF4_dlt);
141  printf("Main: File tof.dat has been read\n");
142  }
143  fclose(F3);
144  UpMat.PlasticThick1/=cos(UpMat.PlasticAngle1*rad);
145  UpMat.PlasticThick2/=cos(UpMat.PlasticAngle2*rad);
146  }
147 
148  /********************* Readout Target parameters:***********************/
149  filePath = gSystem->Getenv("VMCWORKDIR") + TString("/input/target.dat");
150  F1 = fopen(filePath.Data(),"r");
151  if(F1==NULL) printf("Main: File target.dat was not found\n");
152  else
153  {
154  int res;
155  res = fscanf(F1,"%s %s\n",Zeros,UpMat.TarShape);
156  res = fscanf(F1,"%s %lf %lf %lf\n",Zeros,&UpMat.TarXshift,&UpMat.TarYshift,&UpMat.TarZshift);
157  res = fscanf(F1,"%s %s %lf\n",Zeros,UpMat.TarFoilMatter,&UpMat.FoilThick);
158  res = fscanf(F1,"%s %lf\n",Zeros,&UpMat.TarRadius);
159  res = fscanf(F1,"%s %lf\n",Zeros,&UpMat.TarHeight);
160  res = fscanf(F1,"%s %lf\n",Zeros,&UpMat.TarWallThick);
161  res = fscanf(F1,"%s %lf\n",Zeros,&UpMat.TarEntrHoleRad);
162  res = fscanf(F1,"%s %lf\n",Zeros,&UpMat.MeniskSize);
163  res = fscanf(F1,"%s %lf\n",Zeros,&UpMat.TarTemp);
164  res = fscanf(F1,"%s %lf\n",Zeros,&UpMat.TarPress);
165  res = fscanf(F1,"%s %lf\n",Zeros,&UpMat.TarAngle);
166  res = fscanf(F1,"%s %s\n",Zeros,UpMat.HeatScreenAns);
167  res = fscanf(F1,"%s %s %lf\n",Zeros,UpMat.HeatScreenMatter,&UpMat.HeatScreenThick);
168  res = fscanf(F1,"%s %lf\n",Zeros,&UpMat.InHscrRad);
169  res = fscanf(F1,"%s %lf\n",Zeros,&UpMat.HscrWallWidth);
170  res = fscanf(F1,"%s %lf\n",Zeros,&UpMat.HscrHeight);
171  res = fscanf(F1,"%s %lf\n",Zeros,&UpMat.EntrHRad);
172  res = fscanf(F1,"%s %lf\n",Zeros,&UpMat.ExHX);
173  res = fscanf(F1,"%s %lf\n",Zeros,&UpMat.ExHY);
174  printf("Main: File target.dat has been read\n");
175  }
176  fclose(F1);
177  if(UpMat.MeniskSize==0.) UpMat.MeniskSize=0.00001;
178 }
179 //-----------------------------------------------------------------------------
180 void ERHe8EventHeader::ReactionPreparation(){
181  std::cout << " Separate Input and Output channels " << std::endl;
182  char ReaNa[32];
183  char InputChannel[32];
184  char OutputChannel[32];
185  char OutputChannelTemp[32];
186  char Resonance[32];
187  char ResonanceTemp[32];
188  strcpy(ReaNa,ReIN.ReactionName);
189  char* plett;
190  plett = strchr(ReIN.ReactionName,'-');
191  plett+=2;
192  strcpy(OutputChannelTemp,plett);
193  strcpy(ResonanceTemp,plett);
194  plett = strtok(ReIN.ReactionName,"-");
195  strcpy(InputChannel,plett);
196 
197  std::cout << "Define if there is any resonance in the Output channel" << std::endl;
198  int NofPartRes = 0;
199  plett = strchr(OutputChannelTemp,'[');
200  if(plett!=NULL)
201  {
202  plett = strtok(OutputChannelTemp,"[");
203  strcpy(OutputChannel,plett);
204  plett = strtok(NULL,"]");
205  strcat(OutputChannel,plett);
206  plett = strchr(ResonanceTemp,']');plett++;
207  if(plett!=NULL) strcat(OutputChannel,plett);
208  plett = strtok(ResonanceTemp,"[");
209  plett = strtok(NULL,"]");
210  strcpy(ResonanceTemp,plett);
211  plett = strtok(ResonanceTemp,"[");
212  strcpy(ResonanceTemp,plett);
213  plett = strchr(ResonanceTemp,')');
214  if(plett!=NULL)
215  {
216  plett = strtok(ResonanceTemp,")");
217  strcpy(Resonance,plett);
218  plett = strtok(NULL,")");
219  strcat(Resonance,plett);
220  }
221  else strcpy(Resonance,ResonanceTemp);
222  NofPartRes = HowMuchParticles(Resonance);
223  printf("Resonance is %s\n",Resonance);
224  }
225  else strcpy(OutputChannel,OutputChannelTemp);
226 
227  std::cout << "Separate Detected and Unobserved particles in the Output channel" << std::endl;
228  char zero[]="";
229  NofUnObsPart = 0;
230  strcpy(DetectedPart,OutputChannel);
231  plett = strtok(DetectedPart,"()");
232  strcpy(DetectedPart,plett);
233  plett = strchr(OutputChannel,')');
234  plett++;
235  if(strcmp(plett,zero))
236  {
237  plett++;strcpy(UnObservedPart,plett);
238  NofUnObsPart = HowMuchParticles(UnObservedPart);
239  }
240 
241  std::cout << "How much Input and Detected particles" << std::endl;
242  NofInPart = HowMuchParticles(InputChannel);
243  if(NofInPart<2||NofInPart>2) {printf("Wrong number of particles in the Input channel\n");}
244  NofDetPart = HowMuchParticles(DetectedPart);
245  if(NofDetPart==0) {printf("Wrong number of detected particles\n");}
246 
247  std::cout << "Define particles in the input channel" << std::endl;
248 
249  plett = strtok(InputChannel,"+");
250  strcpy(projname,plett);
251  plett = strtok(NULL,"+");
252  strcpy(tarname,plett);
253 }
254 //-----------------------------------------------------------------------------
255 int ERHe8EventHeader::HowMuchParticles(char* str){
256  char xname[32];
257  char* ptr;
258  strcpy(xname,str);
259  int N = 0;
260  do
261  {
262  ptr = strchr(xname,'+');
263  if(ptr!=NULL)
264  {N++;ptr++;strcpy(xname,ptr);}
265  } while(ptr!=NULL);
266  N++;
267  return N;
268 }
269 
270 ClassImp(ERHe8EventHeader)