14 preshStripEnergyCut_(stripEnergyCut), preshSeededNstr_(nStripCut), pathToFiles_(path)
20 string file_pt[5] = {
"20",
"30",
"40",
"50",
"60"};
21 string file_barrel_pt[5] = {
"20",
"30",
"40",
"50",
"60"};
23 string nn_paterns_file =
"";
24 for(
int j=0;
j<Nfiles_EE;
j++) {
25 nn_paterns_file =
"endcapPiZeroDiscriminatorWeights_et"+file_pt[
j]+
".wts";
27 readWeightFile(WFile.fullPath().c_str());
34 for(
int i=0;
i<Indim*Hidden;
i++) I_H_Weight_all.push_back(I_H_Weight[
i]);
35 for(
int i=0;i<Hidden;i++) H_Thresh_all.push_back(H_Thresh[i]);
36 for(
int i=0;i<Outdim*Hidden;i++) H_O_Weight_all.push_back(H_O_Weight[i]);
37 for(
int i=0;i<Outdim;i++) O_Thresh_all.push_back(O_Thresh[i]);
40 for(
int k=0;
k<Nfiles_EB;
k++) {
41 nn_paterns_file =
"barrelPiZeroDiscriminatorWeights_et"+file_barrel_pt[
k]+
".wts";
43 readWeightFile(WFile.fullPath().c_str());
50 for(
int i=0;i<Indim*Hidden;i++) I_H_Weight_all.push_back(I_H_Weight[i]);
51 for(
int i=0;i<Hidden;i++) H_Thresh_all.push_back(H_Thresh[i]);
52 for(
int i=0;i<Outdim*Hidden;i++) H_O_Weight_all.push_back(H_O_Weight[i]);
53 for(
int i=0;i<Outdim;i++) O_Thresh_all.push_back(O_Thresh[i]);
65 vector<float> vout_stripE;
68 if ( rechits_map->size() == 0 ) {
69 edm::LogWarning(
"EndcapPiZeroDiscriminatorAlgo") <<
"RecHitsMap has size 0.";
75 vector<ESDetId> road_2d;
78 int plane = strip.
plane();
80 LogTrace(
"EcalClusters") <<
"EndcapPiZeroDiscriminatorAlgo: findPreshVectors: Preshower Seeded Algorithm - looking for clusters" <<
"n"<<
"findPreshVectors: Preshower is intersected at strip " << strip.
strip() <<
", at plane " << plane ;
84 for(
int i=0;i<11;i++) {
85 vout_stripE.push_back(-100.);
90 road_2d.push_back(strip);
98 LogTrace(
"EcalClusters") <<
"EndcapPiZeroDiscriminatorAlgo:findPreshVectors: Total number of strips in the central road: " << road_2d.size() ;
101 RecHitsMap::iterator final_strip = rechits_map->end();
105 ESDetId last_stripID = final_strip->first;
108 vector<ESDetId>::iterator itID;
109 for (itID = road_2d.begin(); itID != road_2d.end(); itID++) {
110 LogTrace(
"EcalClusters") <<
"EndcapPiZeroDiscriminatorAlgo: findPreshVectors: ID = " << *itID ;
112 RecHitsMap::iterator strip_it = rechits_map->find(*itID);
114 E = strip_it->second.energy();
116 vout_stripE.push_back(E);
117 LogTrace(
"EcalClusters") <<
"EndcapPiZeroDiscriminatorAlgo: findPreshVectors: E = " << E ;
124 vector<float> vout_ElevenStrips_Energy;
125 for(
int i=0;i<11;i++)
127 vout_ElevenStrips_Energy.push_back(0.);
130 for(
unsigned int i=0;i<vout_stripE.size();i++)
132 vout_ElevenStrips_Energy[
i] = vout_stripE.at(i);
136 return vout_ElevenStrips_Energy;
144 RecHitsMap::iterator candidate_tmp = candidate_it;
147 if (candidate_tmp->first == lastID )
148 LogTrace(
"EcalClusters") <<
"EndcapPiZeroDiscriminatorAlgo: goodPi0Strip No such a strip in rechits_map " ;
150 LogTrace(
"EcalClusters") <<
"EndcapPiZeroDiscriminatorAlgo: goodPi0Strip Strip energy " << candidate_it->second.energy() <<
" is below threshold " ;
153 if ( (candidate_tmp->first == lastID ) ||
164 int plane, vector<ESDetId>& vout) {
165 if ( strip ==
ESDetId(0) )
return;
168 LogTrace(
"EcalClusters") <<
"EndcapPiZeroDiscriminatorAlgo: findPi0Road: starts from strip " <<
strip ;
173 LogTrace(
"EcalClusters") <<
"EndcapPiZeroDiscriminatorAlgo: findPi0Road: Go to the East " ;
176 LogTrace(
"EcalClusters") <<
"EndcapPiZeroDiscriminatorAlgo: findPi0Road: East: " << n_east <<
" current strip is " <<
next ;
178 vout.push_back(next);
184 LogTrace(
"EcalClusters") <<
"EndcapPiZeroDiscriminatorAlgo: findPi0Road: Go to the West " ;
188 LogTrace(
"EcalClusters") <<
"findPi0Road: West: " << n_west <<
" current strip is " <<
next ;
190 vout.push_back(next);
194 LogTrace(
"EcalClusters") <<
"EndcapPiZeroDiscriminatorAlgo: findPi0Road: Total number of strips found in the road at 1-st plane is " << n_east+n_west ;
197 else if (plane == 2) {
200 LogTrace(
"EcalClusters") <<
"EndcapPiZeroDiscriminatorAlgo: findPi0Road: Go to the North " ;
203 LogTrace(
"EcalClusters") <<
"EndcapPiZeroDiscriminatorAlgo: findPi0Road: North: " << n_north <<
" current strip is " <<
next;
205 vout.push_back(next);
211 LogTrace(
"EcalClusters") <<
"EndcapPiZeroDiscriminatorAlgo: findPi0Road: Go to the South " ;
215 LogTrace(
"EcalClusters") <<
"EndcapPiZeroDiscriminatorAlgo: findPi0Road: South: " << n_south <<
" current strip is " <<
next ;
217 vout.push_back(next);
221 LogTrace(
"EcalClusters") <<
"EndcapPiZeroDiscriminatorAlgo: findPi0Road: Total number of strips found in the road at 2-nd plane is " << n_south+n_north ;
225 LogTrace(
"EcalClusters") <<
"EndcapPiZeroDiscriminatorAlgo: findPi0Road: Wrong plane number, null cluster will be returned! " ;
244 bool checkinit=
false;
248 weights = fopen(Weights_file,
"r");
249 LogTrace(
"EcalClusters") <<
"EndcapPiZeroDiscriminatorAlgo: I opeded the Weights file = " << Weights_file ;
252 while( !feof(weights) ){
253 fscanf(weights,
"%s", line);
254 if (line[0] ==
'A') {
255 fscanf(weights,
"%d", &
Layers);
256 fscanf(weights,
"%d", &
Indim);
257 fscanf(weights,
"%d", &Hidden);
258 fscanf(weights,
"%d", &
Outdim);
267 }
else if (line[0] ==
'B') {
269 for (
int i = 0; i<
Indim; i++){
274 }
else if (line[0] ==
'C'){
276 for (
int i = 0; i<
Hidden; i++){
277 fscanf(weights,
"%f", &
H_Thresh[i]);
279 }
else if (line[0] ==
'D'){
281 for (
int i = 0; i<Hidden*
Outdim; i++){
284 }
else if (line[0] ==
'E'){
286 for (
int i = 0; i<
Outdim; i++){
287 fscanf(weights,
"%f", &
O_Thresh[i]);
291 else{
edm::LogError(
"EEPi0Discrim")<<
"EndcapPiZeroDiscriminatorAlgo: Not a Net file of Corrupted Net file " << endl;
310 I_SUM =
new float[
Hidden];
314 for(
int k1=0;k1<
Outdim;k1++) OUT[k1]=0.0;
318 for (
int i = 0; i<
Indim; i++){
323 for (
int o1 = 0; o1<
Outdim; o1++) {
328 for (
int o2 = 0; o2<
Outdim; o2++){
332 LogTrace(
"EcalClusters") <<
"EndcapPiZeroDiscriminatorAlgo: getNNoutput :: -> NNout = " << nnout ;
344 return( 1.0 / ( 1.0 +
exp(-2.0*SUM) ) );
358 float pS1_max,
float pS9_max,
float pS25_max,
int EScorr)
362 bool valid_NNinput =
true;
375 for(
int k = 0;
k<11;
k++) {
378 LogTrace(
"EcalClusters") <<
"EndcapPiZeroDiscriminatorAlgo: Oops!!! Preshower Info for strip : " <<
k <<
" of X plane Do not exists" ;
383 LogTrace(
"EcalClusters") <<
"EndcapPiZeroDiscriminatorAlgo: Oops!!! Preshower Info for strip : " <<
k <<
" of Y plane Do not exists" ;
444 for(
int kk1=0;kk1<22;kk1++){
450 float ECAL_norm_factor = 500.;
451 if(pS25_max>500&&pS25_max<=1000) ECAL_norm_factor = 1000;
452 if(pS25_max>1000) ECAL_norm_factor = 7000;
454 input_var[22] = pS1_max/ECAL_norm_factor;
455 input_var[23] = pS9_max/ECAL_norm_factor;
456 input_var[24] = pS25_max/ECAL_norm_factor;
458 LogTrace(
"EcalClusters") <<
"EndcapPiZeroDiscriminatorAlgo: S1/ECAL_norm_factor = " <<
input_var[22];
459 LogTrace(
"EcalClusters") <<
"EndcapPiZeroDiscriminatorAlgo: S9/ECAL_norm_factor = " <<
input_var[23];
460 LogTrace(
"EcalClusters") <<
"EndcapPiZeroDiscriminatorAlgo: S25/ECAL_norm_factor = " <<
input_var[24] ;
464 valid_NNinput =
false;
469 LogTrace(
"EcalClusters") <<
" valid_NNinput = " << valid_NNinput ;
472 return valid_NNinput;
484 double m2,
double cee,
double cep,
double cpp,
485 double s4,
double s6,
double ratio,
486 double xcog,
double ycog)
490 double lam, lam1, lam2;
537 lam=
sqrt((cee -cpp)*(cee -cpp)+4*cep*cep);
538 lam1=(cee + cpp +
lam)/2;
539 lam2=(cee + cpp -
lam)/2;
567 LogTrace(
"EcalClusters") <<
"EndcapPiZeroDiscriminatorAlgo::GetNNoutput :nn_invar_presh = " ;
569 for(
int k1=0;k1<
Indim;k1++) {
577 if(EE_Et<25.0) {sel_wfile = 0;}
578 else if(EE_Et>=25.0 && EE_Et<35.0) {sel_wfile = 1;}
579 else if(EE_Et>=35.0 && EE_Et<45.0) {sel_wfile = 2;}
580 else if(EE_Et>=45.0 && EE_Et<55.0) {sel_wfile = 3;}
581 else {sel_wfile = 4;}
583 LogTrace(
"EcalClusters") <<
"EndcapPiZeroDiscriminatorAlgo: Et_SC = " << EE_Et <<
" and I select Weight file Number = " << sel_wfile ;
588 LogTrace(
"EcalClusters") <<
"EndcapPiZeroDiscriminatorAlgo: ===================> GetNNOutput : NNout = " << nnout ;
608 LogTrace(
"EcalClusters") <<
"EndcapPiZeroDiscriminatorAlgo::GetBarrelNNoutput :nn_invar_presh = " ;
610 for(
int k1=0;k1<
Indim;k1++) {
618 if(EB_Et<25.0) {sel_wfile = 0;}
619 else if(EB_Et>=25.0 && EB_Et<35.0) {sel_wfile = 1;}
620 else if(EB_Et>=35.0 && EB_Et<45.0) {sel_wfile = 2;}
621 else if(EB_Et>=45.0 && EB_Et<55.0) {sel_wfile = 3;}
622 else {sel_wfile = 4;}
623 LogTrace(
"EcalClusters") <<
"EndcapPiZeroDiscriminatorAlgo: E_SC = " << EB_Et <<
" and I select Weight file Number = " << sel_wfile ;
627 LogTrace(
"EcalClusters") <<
"EndcapPiZeroDiscriminatorAlgo: ===================> GetNNOutput : NNout = " << nnout ;
void findPi0Road(ESDetId strip, EcalPreshowerNavigator &theESNav, int plane, std::vector< ESDetId > &vout)
float GetNNOutput(float EE_Et)
void setHome(const T &startingPoint)
set the starting position
void home() const
move the navigator back to the starting point
virtual T west() const
move the navigator west
void calculateBarrelNNInputVariables(float et, double s1, double s9, double s25, double m2, double cee, double cep, double cpp, double s4, double s6, double ratio, double xcog, double ycog)
float Activation_fun(float SUM)
std::vector< float > O_Thresh_all
std::map< DetId, EcalRecHit > RecHitsMap
EndcapPiZeroDiscriminatorAlgo()
virtual T north() const
move the navigator north
double preshStripEnergyCut_
virtual T east() const
move the navigator east
std::vector< float > findPreshVector(ESDetId strip, RecHitsMap *rechits_map, CaloSubdetectorTopology *topology_p)
float getNNoutput(int sel_wfile)
std::vector< float > I_H_Weight_all
float GetBarrelNNOutput(float EB_Et)
bool goodPi0Strip(RecHitsMap::iterator candidate_it, ESDetId lastID)
std::vector< float > H_O_Weight_all
The Signals That Services Can Subscribe To This is based on ActivityRegistry h
Helper function to determine trigger accepts.
void readWeightFile(const char *WFile)
bool calculateNNInputVariables(std::vector< float > &vph1, std::vector< float > &vph2, float pS1_max, float pS9_max, float pS25_max, int EScorr)
virtual T south() const
move the navigator south
std::vector< float > H_Thresh_all