| 43 |
// qtpskip = 0.; |
// qtpskip = 0.; |
| 44 |
lpchi2 = 0.; |
lpchi2 = 0.; |
| 45 |
lpz = 0.; |
lpz = 0.; |
| 46 |
|
lpisotope= 0.; |
| 47 |
lpetot = 0.; |
lpetot = 0.; |
| 48 |
lppskip = 0.; |
lppskip = 0.; |
| 49 |
|
|
| 50 |
memset(calorimetro,0,44*2*sizeof(Float_t)); |
memset(calorimetro,0,44*2*sizeof(Float_t)); |
| 51 |
memset(spessore,0,4*sizeof(Float_t)); |
memset(spessore,0,4*sizeof(Float_t)); |
| 52 |
memset(estremi,0,2*2*sizeof(Float_t)); |
memset(estremi,0,2*2*sizeof(Float_t)); |
| 73 |
// printf(" plane not used for truncated mean %f: \n", qtpskip); |
// printf(" plane not used for truncated mean %f: \n", qtpskip); |
| 74 |
printf(" chi 2 from loop %f: \n", lpchi2); |
printf(" chi 2 from loop %f: \n", lpchi2); |
| 75 |
printf(" Z from loop %f: \n", lpz); |
printf(" Z from loop %f: \n", lpz); |
| 76 |
|
printf(" isotope from loop %f: \n", lpisotope); |
| 77 |
printf(" energy from loop %f: \n", lpetot); |
printf(" energy from loop %f: \n", lpetot); |
| 78 |
printf(" plane not used for loop %f: \n", lppskip); |
printf(" plane not used for loop %f: \n", lppskip); |
| 79 |
printf(" ndf: %i \n",ndf); |
printf(" ndf: %i \n",ndf); |
| 92 |
}; |
}; |
| 93 |
|
|
| 94 |
|
|
| 95 |
void CaloBragg::CleanPlanes(Float_t epiano[22][2]){ |
void CaloBragg::CleanPlanes(Float_t epiano[22][2], Bool_t zpiano[22][2]){ |
| 96 |
// return; |
// return; |
| 97 |
Int_t hitplanes = 0; |
Int_t hitplanes = 0; |
| 98 |
|
Float_t f5 = 0.; |
| 99 |
for (Int_t i = 0; i<22; i++){ |
for (Int_t i = 0; i<22; i++){ |
| 100 |
for (Int_t j = 1; j>=0; j--){ |
for (Int_t j = 1; j>=0; j--){ |
| 101 |
if ( epiano[i][j] > 0.7 ) hitplanes++; |
zpiano[i][j] = false; |
| 102 |
|
if ( epiano[i][j] > 0.7 ){ |
| 103 |
|
if ( hitplanes < 100 ) f5 += epiano[i][j]; |
| 104 |
|
hitplanes++; |
| 105 |
|
}; |
| 106 |
}; |
}; |
| 107 |
}; |
}; |
| 108 |
|
Int_t atl5 = TMath::Min(hitplanes,100); |
| 109 |
|
atl5 = TMath::Max(atl5,1); |
| 110 |
Float_t lowlim = 0.85; |
Float_t lowlim = 0.85; |
| 111 |
|
//Float_t lowlim = 1.; |
| 112 |
Float_t dedxone = 0.; |
Float_t dedxone = 0.; |
| 113 |
Float_t step1 = 0.8*L2->GetCaloLevel2()->qtot/(Float_t)hitplanes; |
// Float_t step1 = 0.8*L2->GetCaloLevel2()->qtot/(Float_t)hitplanes; |
| 114 |
while ( dedxone < step1 ){ |
Float_t step1 = 0.8*f5/atl5; |
| 115 |
|
// while ( dedxone < step1 ){ |
| 116 |
for (Int_t i = 0; i<22; i++){ |
for (Int_t i = 0; i<22; i++){ |
| 117 |
for (Int_t j = 1; j>=0; j--){ |
for (Int_t j = 1; j>=0; j--){ |
| 118 |
|
if (debug) printf("Acleanplanes: i %i j %i step1 %f dedxone %f epiano[i][j] %f \n",i,j,step1,dedxone,epiano[i][j]); |
| 119 |
if ( epiano[i][j] >= step1 && dedxone < 0.7 ) dedxone = epiano[i][j]; |
if ( epiano[i][j] >= step1 && dedxone < 0.7 ) dedxone = epiano[i][j]; |
| 120 |
|
if ( dedxone >= step1 ) break; // new |
| 121 |
}; |
}; |
| 122 |
|
if ( dedxone >= step1 ) break; // new |
| 123 |
}; |
}; |
| 124 |
} |
// } |
| 125 |
if ( dedxone < 0.7 ){ |
if ( dedxone < 0.7 ){ // here we could have instead while dedxone == 0. ... perhaps better... |
| 126 |
for (Int_t i = 0; i<22; i++){ |
for (Int_t i = 0; i<22; i++){ |
| 127 |
for (Int_t j = 1; j>=0; j--){ |
for (Int_t j = 1; j>=0; j--){ |
| 128 |
|
if (debug) printf("Bcleanplanes dedxone < 0.7: i %i j %i step1 %f dedxone %f epiano[i][j] %f \n",i,j,step1,dedxone,epiano[i][j]); |
| 129 |
if ( epiano[i][j] > 0. && dedxone < 0.7 ) dedxone = epiano[i][j]; |
if ( epiano[i][j] > 0. && dedxone < 0.7 ) dedxone = epiano[i][j]; |
| 130 |
|
if ( dedxone >= 0.7 ) break; // new |
| 131 |
}; |
}; |
| 132 |
|
if ( dedxone >= 0.7 ) break; // new |
| 133 |
}; |
}; |
| 134 |
} |
} |
| 135 |
// |
// |
| 140 |
for (Int_t i = 0; i<22; i++){ |
for (Int_t i = 0; i<22; i++){ |
| 141 |
for (Int_t j = 1; j>=0; j--){ |
for (Int_t j = 1; j>=0; j--){ |
| 142 |
if ( epiano[i][j] < dedxone*lowlim ){ |
if ( epiano[i][j] < dedxone*lowlim ){ |
| 143 |
// printf(" %i %i epiano %f limit %f nulliferus %i nullius %i \n",i,j,epiano[i][j],dedxone*lowlim,nulliferus,nullius); |
if ( debug ) printf("Ccleanplanes: %i %i epiano %f limit %f nulliferus %i nullius %i \n",i,j,epiano[i][j],dedxone*lowlim,nulliferus,nullius); |
| 144 |
epiano[i][j] = 0.; |
// epiano[i][j] = 0.; |
| 145 |
|
zpiano[i][j] = true; |
| 146 |
|
if ( epiano[i][j] < dedxone*0.05 ) epiano[i][j] = 0.; |
| 147 |
} else { |
} else { |
| 148 |
//x printf(" %i %i epiano %f limit %f nulliferus %i nullius %i \n",i,j,epiano[i][j],dedxone*lowlim,nulliferus,nullius); |
if ( debug ) printf("Dcleanplanes else: %i %i epiano %f limit %f nulliferus %i nullius %i \n",i,j,epiano[i][j],dedxone*lowlim,nulliferus,nullius); |
| 149 |
nulliferus = 0; |
nulliferus = 0; |
| 150 |
revulsera = true; |
revulsera = true; |
| 151 |
}; |
}; |
| 152 |
if ( epiano[i][j] < 0.7 && revulsera ) nulliferus++; |
// if ( epiano[i][j] < 0.7 && revulsera ) nulliferus++; |
| 153 |
|
if ( (zpiano[i][j] || epiano[i][j] < 0.7 ) && revulsera ) nulliferus++; |
| 154 |
if ( nulliferus > 10 ) nullius = true; |
if ( nulliferus > 10 ) nullius = true; |
| 155 |
if ( nullius ) epiano[i][j] = 0.; |
// if ( nullius ) epiano[i][j] = 0.; |
| 156 |
|
if ( nullius ) zpiano[i][j] = true; |
| 157 |
}; |
}; |
| 158 |
}; |
}; |
| 159 |
|
|
| 210 |
// |
// |
| 211 |
}; |
}; |
| 212 |
// |
// |
| 213 |
this->CleanPlanes(*&epiano); |
Bool_t zpiano[22][2]; |
| 214 |
|
this->CleanPlanes(*&epiano, *&zpiano); |
| 215 |
// |
// |
| 216 |
PamTrack *ptrack = 0; |
PamTrack *ptrack = 0; |
| 217 |
CaloTrkVar *track = 0; |
CaloTrkVar *track = 0; |
| 240 |
if(L2->GetCaloLevel2()->npcfit[0]==0 && L2->GetCaloLevel2()->npcfit[1]==0 && L2->GetCaloLevel2()->npcfit[2]==0 && L2->GetCaloLevel2()->npcfit[3]==0) return;// controllo sulla traccia nel calorimetro |
if(L2->GetCaloLevel2()->npcfit[0]==0 && L2->GetCaloLevel2()->npcfit[1]==0 && L2->GetCaloLevel2()->npcfit[2]==0 && L2->GetCaloLevel2()->npcfit[3]==0) return;// controllo sulla traccia nel calorimetro |
| 241 |
|
|
| 242 |
// |
// |
| 243 |
|
Bool_t zcalo[44]; |
| 244 |
for(Int_t p=0; p<22; p++){ |
for(Int_t p=0; p<22; p++){ |
| 245 |
for(Int_t v=0; v<2; v++){ |
for(Int_t v=0; v<2; v++){ |
| 246 |
/*per usare traccia non del calo camboare cibar*/ |
/*per usare traccia non del calo camboare cibar*/ |
| 247 |
calorimetro[(2*p)+1-v][0] = L2->GetCaloLevel2()->cibar[p][v];//strip attraversata |
calorimetro[(2*p)+1-v][0] = L2->GetCaloLevel2()->cibar[p][v];//strip attraversata |
| 248 |
calorimetro[(2*p)+1-v][1] = epiano[p][v]; //energia del piano //(epiano[p][v])/0.89 |
calorimetro[(2*p)+1-v][1] = epiano[p][v]; //energia del piano //(epiano[p][v])/0.89 |
| 249 |
|
zcalo[(2*p)+1-v] = zpiano[p][v]; |
| 250 |
|
if ( debug ) printf(" idx %i %f %i \n",(2*p)+1-v,epiano[p][v], zpiano[p][v]); |
| 251 |
}; |
}; |
| 252 |
}; |
}; |
| 253 |
|
|
| 311 |
// |
// |
| 312 |
|
|
| 313 |
Float_t lastok = 0.; |
Float_t lastok = 0.; |
| 314 |
|
// if ( false ){ |
| 315 |
// Bool_t goback = false; |
// Bool_t goback = false; |
| 316 |
for ( int o = 0; o < estremi[1][0]; o++ ){ |
for ( int o = 0; o < estremi[1][0]; o++ ){ |
| 317 |
// |
// |
| 318 |
if ( calorimetro[o][1] > 0.7 ) lastok = calorimetro[o][1]; |
if (debug) printf(" goforth1: o %i calo %f lastok %f \n",o,calorimetro[o][1],lastok); |
| 319 |
if ( calorimetro[o][1] < 0.7 && lastok > 0. ) calorimetro[o][1] = lastok; |
if ( calorimetro[o][1] > 0.7 && !zcalo[o] ) lastok = calorimetro[o][1]; |
| 320 |
|
if ( (zcalo[o] || calorimetro[o][1] < 0.7) && lastok > 0. ){ |
| 321 |
|
if ( fabs(calorimetro[o][1]-lastok)/calorimetro[o][1] > 0.5 ) { |
| 322 |
|
if (debug) printf(" goforthchange %f %f \n",calorimetro[o][1],lastok); |
| 323 |
|
calorimetro[o][1] = lastok; |
| 324 |
|
if (debug) printf(" goforthchang+ %f %f \n",calorimetro[o][1],lastok); |
| 325 |
|
} |
| 326 |
|
} |
| 327 |
|
if (debug) printf(" goforth2: o %i calo %f lastok %f \n",o,calorimetro[o][1],lastok); |
| 328 |
// if ( calorimetro[o][1] < 0.7 ) goback = true; |
// if ( calorimetro[o][1] < 0.7 ) goback = true; |
| 329 |
// |
// |
| 330 |
}; |
}; |
| 332 |
// if ( goback ){ |
// if ( goback ){ |
| 333 |
for ( int o = estremi[1][0]; o >= 0; o-- ){ |
for ( int o = estremi[1][0]; o >= 0; o-- ){ |
| 334 |
// |
// |
| 335 |
// printf(" goback1: o %i calo %f lastok %f \n",o,calorimetro[o][1],lastok); |
if (debug) printf(" goback1: o %i calo %f lastok %f \n",o,calorimetro[o][1],lastok); |
| 336 |
if ( o < estremi[1][0] && calorimetro[o][1] > calorimetro[o+1][1]*1.2 && lastok > 0. ) calorimetro[o][1] = lastok; |
if ( o < estremi[1][0] && calorimetro[o][1] > calorimetro[o+1][1]*1.2 && lastok > 0. ) calorimetro[o][1] = lastok; |
| 337 |
if ( calorimetro[o][1] > 0.7 ) lastok = calorimetro[o][1]; |
if ( calorimetro[o][1] > 0.7 && !zcalo[o] ) lastok = calorimetro[o][1]; |
| 338 |
if ( calorimetro[o][1] < 0.7 && lastok > 0. ) calorimetro[o][1] = lastok; |
if ( (zcalo[o] || calorimetro[o][1] < 0.7) && lastok > 0. ){ |
| 339 |
// printf(" goback2: o %i calo %f lastok %f \n",o,calorimetro[o][1],lastok); |
if ( fabs(calorimetro[o][1]-lastok)/calorimetro[o][1] > 0.5 ) { |
| 340 |
|
if ( debug ) printf(" gobackchange \n"); |
| 341 |
|
calorimetro[o][1] = lastok; |
| 342 |
|
} |
| 343 |
|
} |
| 344 |
|
if (debug) printf(" goback2: o %i calo %f lastok %f \n",o,calorimetro[o][1],lastok); |
| 345 |
// |
// |
| 346 |
}; |
}; |
| 347 |
// }; |
// }; |
| 348 |
|
//} |
| 349 |
|
|
| 350 |
if ( startZero ) { |
if ( startZero ) { |
| 351 |
estremi[0][0] = 0.; |
estremi[0][0] = 0.; |
| 384 |
|
|
| 385 |
Float_t dEpianiloop[44]; |
Float_t dEpianiloop[44]; |
| 386 |
Int_t tz1=(Int_t)lpz; |
Int_t tz1=(Int_t)lpz; |
| 387 |
Enetrack(&tz1, &lpetot, &estremi[0][0],&estremi[1][0], dEpianiloop);//calcola rilascio energetico sui piani da loop |
Int_t ti1=(Int_t)lpisotope; |
| 388 |
|
|
| 389 |
|
Enetrack(&tz1, &ti1 , &lpetot, &estremi[0][0],&estremi[1][0], dEpianiloop);//calcola rilascio energetico sui piani da loop |
| 390 |
|
|
| 391 |
|
|
| 392 |
Float_t integ = 0.; |
Float_t integ = 0.; |
| 415 |
|
|
| 416 |
Process(); |
Process(); |
| 417 |
|
|
| 418 |
this->Draw(0.,0.); |
this->Draw(0.,0.,0.); |
| 419 |
|
|
| 420 |
} |
} |
| 421 |
|
|
| 422 |
void CaloBragg::Draw(Int_t Z, Float_t enetot){ |
void CaloBragg::Draw(Int_t Z, Int_t isotope, Float_t enetot){ |
| 423 |
|
|
| 424 |
// Float_t dEpianimean[44]; |
// Float_t dEpianimean[44]; |
| 425 |
Float_t dEpianiloop[44]; |
Float_t dEpianiloop[44]; |
| 426 |
Float_t Depth[44]; |
Float_t Depth[44]; |
| 427 |
// Int_t tz=(Int_t)qtz; |
// Int_t tz=(Int_t)qtz; |
| 428 |
Int_t tz1= Z; |
Int_t tz1= Z; |
| 429 |
|
Int_t ti1= isotope; |
| 430 |
Float_t enet = enetot; |
Float_t enet = enetot; |
| 431 |
// Float_t enet = lpetot; |
// Float_t enet = lpetot; |
| 432 |
|
|
| 458 |
|
|
| 459 |
} else { |
} else { |
| 460 |
tz1=(Int_t)lpz; |
tz1=(Int_t)lpz; |
| 461 |
|
ti1=(Int_t)lpisotope; |
| 462 |
enet = lpetot; |
enet = lpetot; |
| 463 |
// Enetrack(&tz, &qtetot, &estremi[0][0],&estremi[1][0], dEpianimean);//calcola rilascio energetico sui piani da media troncata |
// Enetrack(&tz, &qtetot, &estremi[0][0],&estremi[1][0], dEpianimean);//calcola rilascio energetico sui piani da media troncata |
| 464 |
|
|
| 465 |
} |
} |
| 466 |
Enetrack(&tz1, &enet, &estremi[0][0],&estremi[1][0], dEpianiloop);//calcola rilascio energetico sui piani da loop |
Enetrack(&tz1, &ti1, &enet, &estremi[0][0],&estremi[1][0], dEpianiloop);//calcola rilascio energetico sui piani da loop |
| 467 |
|
|
| 468 |
Float_t sp= spessore[0]*spessore[1]; |
Float_t sp= spessore[0]*spessore[1]; |
| 469 |
for(Int_t i=0;i<44;i++)Depth[i]=i*sp; |
for(Int_t i=0;i<44;i++)Depth[i]=i*sp; |
| 494 |
TString thid2 = Form("hCaloBragg2"); |
TString thid2 = Form("hCaloBragg2"); |
| 495 |
TH2F *th2 = dynamic_cast<TH2F*>(gDirectory->FindObject(thid2)); |
TH2F *th2 = dynamic_cast<TH2F*>(gDirectory->FindObject(thid2)); |
| 496 |
if ( th2 ) th2->Delete(); |
if ( th2 ) th2->Delete(); |
| 497 |
th2 = new TH2F(thid2,thid2,300,-0.5,300.,1000,0.,150.); |
th2 = new TH2F(thid2,thid2,300,-0.5,300.,1000,0.,25.); //150 |
| 498 |
th2->SetMarkerStyle(20); |
th2->SetMarkerStyle(20); |
| 499 |
th2->SetMarkerColor(kRed); |
th2->SetMarkerColor(kRed); |
| 500 |
// |
// |
| 501 |
TString thid3 = Form("hCaloBragg3"); |
TString thid3 = Form("hCaloBragg3"); |
| 502 |
TH2F *th3 = dynamic_cast<TH2F*>(gDirectory->FindObject(thid3)); |
TH2F *th3 = dynamic_cast<TH2F*>(gDirectory->FindObject(thid3)); |
| 503 |
if ( th3 ) th3->Delete(); |
if ( th3 ) th3->Delete(); |
| 504 |
th3 = new TH2F(thid3,thid3,300,-0.5,300.,1000,0.,150.); |
th3 = new TH2F(thid3,thid3,300,-0.5,300.,1000,0.,25.);//150. |
| 505 |
th3->SetMarkerStyle(20); |
th3->SetMarkerStyle(20); |
| 506 |
th3->SetMarkerColor(kBlue); |
th3->SetMarkerColor(kBlue); |
| 507 |
|
|
| 536 |
|
|
| 537 |
void CaloBragg::LoadParam(){ |
void CaloBragg::LoadParam(){ |
| 538 |
|
|
| 539 |
// |
// elem[Z-1][isotop] 0 is the most common one |
| 540 |
elem[0] = 1.00794; //H 1 |
// |
| 541 |
elem[1] = 4.002602; //He 2 |
|
| 542 |
elem[2] = 6.941; //Li 3 |
elem[0][0] = 1.00782; //H 1 |
| 543 |
elem[3] = 9.012182;//Be 4 |
elem[0][1] = 2.01410; // 2H (Isotope) |
| 544 |
elem[4] = 10.811; //B 5 |
elem[0][2] = -1.; |
| 545 |
elem[5] = 12.0107; //C 6 |
elem[0][3] = -1.; |
| 546 |
elem[6] = 14.00674;//N 7 |
elem[0][4] = -1.; |
| 547 |
elem[7] = 15.9994; //O 8 |
elem[0][5] = -1.; |
| 548 |
elem[8] = 18.9984032; //F 9 |
elem[0][6] = -1.; |
| 549 |
elem[9] = 20.1797; //Ne 10 |
|
| 550 |
elem[10] = 22.98977;//Na 11 |
elem[1][0] = 4.002602; //He 2 |
| 551 |
elem[11] = 24.3050; //Mg 12 |
elem[1][1] = 3.016029; // 3He (Isotope) |
| 552 |
elem[12] = 26.9815; //Al 13 |
elem[1][2] = -1.; |
| 553 |
elem[13] = 28.0855; //Si 14 |
elem[1][3] = -1.; |
| 554 |
elem[14] = 30.974; //P 15 |
elem[1][4] = -1.; |
| 555 |
elem[15] = 32.066; //S 16 |
elem[1][5] = -1.; |
| 556 |
elem[16] = 35.4527; //Cl 17 |
elem[1][6] = -1.; |
| 557 |
elem[17] = 39.948; //Ar 18 |
|
| 558 |
elem[18] = 39.0983; //K 19 |
elem[2][0] = 7.016004; //Li 3 |
| 559 |
elem[19] = 40.078; //Ca 20 |
elem[2][1] = 6.015123; //6Li (Isotope) |
| 560 |
elem[20] = 44.95591;//Sc 21 |
elem[2][2] = -1.; |
| 561 |
elem[21] = 47.867; //Ti 22 |
elem[2][3] = -1.; |
| 562 |
elem[22] = 50.9415; //V 23 |
elem[2][4] = -1.; |
| 563 |
elem[23] = 51.9961; //Cr 24 |
elem[2][5] = -1.; |
| 564 |
elem[24] = 54.938049;//Mn 25 |
elem[2][6] = -1.; |
| 565 |
elem[25] = 55.845; //Fe 26 |
|
| 566 |
elem[26] = 58.9332; //Co 27 |
elem[3][0] = 9.012182; //Be 4 |
| 567 |
elem[27] = 58.6934; //Ni 28 |
elem[3][1] = 10.01353; //10Be (Isotope) (most stable) |
| 568 |
elem[28] = 63.546; //Cu 29 |
elem[3][2] = 7.01693; //9Be no EC in space? |
| 569 |
elem[29] = 65.39; //Zn 30 |
elem[3][3] = -1.; |
| 570 |
elem[30] = 69.723; //Ga 31 |
elem[3][4] = -1.; |
| 571 |
elem[31] = 72.61; //Ge 32 |
elem[3][5] = -1.; |
| 572 |
|
elem[3][6] = -1.; |
| 573 |
|
|
| 574 |
|
elem[4][0] = 11.00930; //B 5 |
| 575 |
|
elem[4][1] = 10.01294; //10B (Isotope) |
| 576 |
|
elem[4][2] = -1.; |
| 577 |
|
elem[4][3] = -1.; |
| 578 |
|
elem[4][4] = -1.; |
| 579 |
|
elem[4][5] = -1.; |
| 580 |
|
elem[4][5] = -1.; |
| 581 |
|
|
| 582 |
|
elem[5][0] = 12.0107; //C 6 |
| 583 |
|
elem[5][1] = 13.00335; //13C (Isotope) |
| 584 |
|
elem[5][2] = -1.; |
| 585 |
|
elem[5][3] = -1.; |
| 586 |
|
elem[5][4] = -1.; |
| 587 |
|
elem[5][5] = -1.; |
| 588 |
|
elem[5][5] = -1.; |
| 589 |
|
|
| 590 |
|
elem[6][0] = 14.00674; //N 7 |
| 591 |
|
elem[6][1] = 15.00011; //15N (Isotope) |
| 592 |
|
elem[6][2] = -1.; |
| 593 |
|
elem[6][3] = -1.; |
| 594 |
|
elem[6][4] = -1.; |
| 595 |
|
elem[6][5] = -1.; |
| 596 |
|
elem[6][5] = -1.; |
| 597 |
|
|
| 598 |
|
elem[7][0] = 15.99491; //O 8 |
| 599 |
|
elem[7][1] = 17.99916; //18O (Isotope) |
| 600 |
|
elem[7][2] = 16.99916; //17O (Isotope) |
| 601 |
|
elem[7][3] = -1.; |
| 602 |
|
elem[7][4] = -1.; |
| 603 |
|
elem[7][5] = -1.; |
| 604 |
|
elem[7][5] = -1.; |
| 605 |
|
|
| 606 |
|
elem[8][0] = 18.99840; //F 9 |
| 607 |
|
elem[8][1] = -1.; |
| 608 |
|
elem[8][2] = -1.; |
| 609 |
|
elem[8][3] = -1.; |
| 610 |
|
elem[8][4] = -1.; |
| 611 |
|
elem[8][5] = -1.; |
| 612 |
|
elem[8][5] = -1.; |
| 613 |
|
|
| 614 |
|
elem[9][0] = 19.99244; //Ne 10 |
| 615 |
|
elem[9][1] = 21.99138; //22Ne (Isotope) |
| 616 |
|
elem[9][2] = 20.99384; //21Ne 10 |
| 617 |
|
elem[9][3] = -1.; |
| 618 |
|
elem[9][4] = -1.; |
| 619 |
|
elem[9][5] = -1.; |
| 620 |
|
elem[9][6] = -1.; |
| 621 |
|
|
| 622 |
|
elem[10][0] = 22.98977; //Na 11 |
| 623 |
|
elem[10][1] = 21.99444; //22Na (Isotope) (most stable) |
| 624 |
|
elem[10][2] = -1.; |
| 625 |
|
elem[10][3] = -1.; |
| 626 |
|
elem[10][4] = -1.; |
| 627 |
|
elem[10][5] = -1.; |
| 628 |
|
elem[10][6] = -1.; |
| 629 |
|
|
| 630 |
|
elem[11][0] = 23.98504; //Mg 12 |
| 631 |
|
elem[11][1] = 25.98259; //26Mg (Isotope) |
| 632 |
|
elem[11][2] = 24.98504; //25Mg (Isotope) |
| 633 |
|
elem[11][3] = -1.; |
| 634 |
|
elem[11][4] = -1.; |
| 635 |
|
elem[11][5] = -1.; |
| 636 |
|
elem[11][6] = -1.; |
| 637 |
|
|
| 638 |
|
elem[12][0] = 26.98154; //Al 13 |
| 639 |
|
elem[12][1] = 25.98489; //26Al (Isotope) (most stable) |
| 640 |
|
elem[12][2] = -1.; |
| 641 |
|
elem[12][3] = -1.; |
| 642 |
|
elem[12][4] = -1.; |
| 643 |
|
elem[12][5] = -1.; |
| 644 |
|
elem[12][6] = -1.; |
| 645 |
|
|
| 646 |
|
elem[13][0] = 27.97692; //Si 14 |
| 647 |
|
elem[13][1] = 28.97649; //29Si (Isotope) |
| 648 |
|
elem[13][2] = 29.97377; //30Si (Isotope) |
| 649 |
|
elem[13][3] = -1.; |
| 650 |
|
elem[13][4] = -1.; |
| 651 |
|
elem[13][5] = -1.; |
| 652 |
|
elem[13][6] = -1.; |
| 653 |
|
|
| 654 |
|
elem[14][0] = 30.97376; //P 15 |
| 655 |
|
elem[14][1] = -1.; |
| 656 |
|
elem[14][2] = -1.; |
| 657 |
|
elem[14][3] = -1.; |
| 658 |
|
elem[14][4] = -1.; |
| 659 |
|
elem[14][5] = -1.; |
| 660 |
|
elem[14][6] = -1.; |
| 661 |
|
|
| 662 |
|
elem[15][0] = 31.97207; //S 16 |
| 663 |
|
elem[15][1] = 33.96787; //34S (Isotope) |
| 664 |
|
elem[15][2] = 32.97146; //33S (Isotope) |
| 665 |
|
elem[15][3] = 35.96708; //36S (Isotope) |
| 666 |
|
elem[15][4] = -1.; |
| 667 |
|
elem[15][5] = -1.; |
| 668 |
|
elem[15][6] = -1.; |
| 669 |
|
|
| 670 |
|
elem[16][0] = 34.96885; //Cl 17 |
| 671 |
|
elem[16][1] = 36.96831; //37Cl 17 |
| 672 |
|
elem[16][2] = 35.96890; //36Cl (Isotope) |
| 673 |
|
elem[16][3] = -1.; |
| 674 |
|
elem[16][4] = -1.; |
| 675 |
|
elem[16][5] = -1.; |
| 676 |
|
elem[16][6] = -1.; |
| 677 |
|
|
| 678 |
|
elem[17][0] = 39.962383; //Ar 18 |
| 679 |
|
elem[17][1] = 35.967545; //36Ar (Isotope) |
| 680 |
|
elem[17][2] = 37.962732; //38Ar (Isotope) |
| 681 |
|
elem[17][3] = 38.964313; //39Ar (Isotope) |
| 682 |
|
elem[17][4] = -1.; |
| 683 |
|
elem[17][5] = -1.; |
| 684 |
|
elem[17][6] = -1.; |
| 685 |
|
|
| 686 |
|
elem[18][0] = 38.963707; //K 19 |
| 687 |
|
elem[18][1] = 40.961825; //41K (Isotope) |
| 688 |
|
elem[18][2] = 39.963998; //40K (Isotope) |
| 689 |
|
elem[18][3] = -1.; |
| 690 |
|
elem[18][4] = -1.; |
| 691 |
|
elem[18][5] = -1.; |
| 692 |
|
elem[18][6] = -1.; |
| 693 |
|
|
| 694 |
|
elem[19][0] = 39.962590; //Ca 20 |
| 695 |
|
elem[19][1] = 43.955482; //44Ca (Isotope) |
| 696 |
|
elem[19][2] = 41.958618; //42Ca (Isotope) |
| 697 |
|
elem[19][3] = 42.958767; //43Ca (Isotope) |
| 698 |
|
elem[19][4] = 45.953693; //46Ca (Isotope) |
| 699 |
|
elem[19][5] = 40.962278; //41Ca (Isotope) |
| 700 |
|
elem[19][6] = -1.; |
| 701 |
|
|
| 702 |
|
elem[20][0] = 44.955912;//Sc 21 |
| 703 |
|
elem[20][1] = -1.; |
| 704 |
|
elem[20][2] = -1.; |
| 705 |
|
elem[20][3] = -1.; |
| 706 |
|
elem[20][4] = -1.; |
| 707 |
|
elem[20][5] = -1.; |
| 708 |
|
elem[20][6] = -1.; |
| 709 |
|
|
| 710 |
|
elem[21][0] = 47.947946; //Ti 22 |
| 711 |
|
elem[21][1] = 45.952632; //46Ti (Isotope) |
| 712 |
|
elem[21][2] = 46.951763; //47Ti (Isotope) |
| 713 |
|
elem[21][3] = 48.947870; //49Ti (Isotope) |
| 714 |
|
elem[21][4] = 49.944791; //50Ti (Isotope) |
| 715 |
|
elem[21][5] = 43.959690; //44Ti (Isotope) (half life 60y) |
| 716 |
|
elem[21][6] = -1.; |
| 717 |
|
|
| 718 |
|
elem[22][0] = 50.943960; //V 23 |
| 719 |
|
elem[22][1] = 49.947158; //50V (Isotope) |
| 720 |
|
elem[22][2] = -1.; |
| 721 |
|
elem[22][3] = -1.; |
| 722 |
|
elem[22][4] = -1.; |
| 723 |
|
elem[22][5] = -1.; |
| 724 |
|
elem[22][6] = -1.; |
| 725 |
|
|
| 726 |
|
elem[23][0] = 51.940507; //Cr 24 |
| 727 |
|
elem[23][1] = 52.940649; //53Cr (Isotope) |
| 728 |
|
elem[23][2] = 49.946044; //50Cr (Isotope) |
| 729 |
|
elem[23][3] = 53.938880; //54Cr (Isotope) |
| 730 |
|
elem[23][4] = -1.; |
| 731 |
|
elem[23][5] = -1.; |
| 732 |
|
elem[23][6] = -1.; |
| 733 |
|
|
| 734 |
|
elem[24][0] = 54.938049;//Mn 25 |
| 735 |
|
elem[24][1] = 52.941290;//53Mn (Isotope) |
| 736 |
|
elem[24][2] = -1.; |
| 737 |
|
elem[24][3] = -1.; |
| 738 |
|
elem[24][4] = -1.; |
| 739 |
|
elem[24][5] = -1.; |
| 740 |
|
elem[24][6] = -1.; |
| 741 |
|
|
| 742 |
|
elem[25][0] = 55.934937; //Fe 26 |
| 743 |
|
elem[25][1] = 53.939610; //54Fe (Isotope) |
| 744 |
|
elem[25][2] = 56.935394; //57Fe (Isotope) |
| 745 |
|
elem[25][3] = 57.933276; //58Fe (Isotope) |
| 746 |
|
elem[25][4] = 59.934072; //58Fe (Isotope) |
| 747 |
|
|
| 748 |
|
elem[26][0] = 58.933195; //Co 27 |
| 749 |
|
elem[26][1] = 59.933817; //60Co (Isotope) |
| 750 |
|
elem[26][2] = -1.; |
| 751 |
|
elem[26][3] = -1.; |
| 752 |
|
elem[26][4] = -1.; |
| 753 |
|
elem[26][5] = -1.; |
| 754 |
|
elem[26][6] = -1.; |
| 755 |
|
|
| 756 |
|
|
| 757 |
|
elem[27][0] = 57.935343; //Ni 28 |
| 758 |
|
elem[27][1] = 61.928345; //62Ni (Isotope) |
| 759 |
|
elem[27][2] = 59.930786; //60Ni (Isotope) |
| 760 |
|
elem[27][3] = 60.931056; //61Ni (Isotope) |
| 761 |
|
elem[27][4] = 63.927966; //64Ni (Isotope) |
| 762 |
|
elem[27][5] = 58.934346; //59Ni (Isotope) |
| 763 |
|
elem[27][6] = -1.; |
| 764 |
|
|
| 765 |
|
elem[28][0] = 62.929597; //Cu 29 |
| 766 |
|
elem[28][1] = 64.927789; //65Cu (Isotope) |
| 767 |
|
elem[28][2] = -1.; |
| 768 |
|
elem[28][3] = -1.; |
| 769 |
|
elem[28][4] = -1.; |
| 770 |
|
elem[28][5] = -1.; |
| 771 |
|
elem[28][6] = -1.; |
| 772 |
|
|
| 773 |
|
elem[29][0] = 63.929142; //Zn 30 |
| 774 |
|
elem[29][1] = 65.926033; //66Zn (Isotope) |
| 775 |
|
elem[29][2] = 67.924844; //68Zn (Isotope) |
| 776 |
|
elem[29][3] = 66.927127; //67Zn (Isotope) |
| 777 |
|
elem[29][4] = 69.925319; //70Zn (Isotope) |
| 778 |
|
elem[29][5] = -1.; |
| 779 |
|
elem[29][6] = -1.; |
| 780 |
|
|
| 781 |
|
elem[30][0] = 68.925573; //Ga 31 |
| 782 |
|
elem[30][1] = 70.924701; //71Ga (Isotope) |
| 783 |
|
elem[30][2] = -1.; |
| 784 |
|
elem[30][3] = -1.; |
| 785 |
|
elem[30][4] = -1.; |
| 786 |
|
elem[30][5] = -1.; |
| 787 |
|
elem[30][6] = -1.; |
| 788 |
|
|
| 789 |
|
elem[31][0] = 73.921177; //Ge 32 |
| 790 |
|
elem[31][1] = 71.922075; //72Ge (Isotope) |
| 791 |
|
elem[31][2] = 69.924247; //70Ge (Isotope) |
| 792 |
|
elem[31][3] = 75.921403; //76Ge (Isotope) |
| 793 |
|
elem[31][4] = 73.923459; //73Ge (Isotope) |
| 794 |
|
elem[31][5] = -1.; |
| 795 |
|
elem[31][6] = -1.; |
| 796 |
|
|
| 797 |
|
|
| 798 |
//parametri calorimetro |
//parametri calorimetro |
| 919 |
|
|
| 920 |
|
|
| 921 |
|
|
| 922 |
void CaloBragg::ELOSS(Float_t *dx, Int_t *Z, Float_t *Etot, Float_t *out, Float_t II){ |
void CaloBragg::ELOSS(Float_t *dx, Int_t *Z, Int_t *isotope, Float_t *Etot, Float_t *out, Float_t II){ |
| 923 |
|
|
| 924 |
/*perdita di energia per ioni pesanti (come da routine geant)*/ |
/*perdita di energia per ioni pesanti (come da routine geant)*/ |
| 925 |
// in : dx => spessore g/cm2 |
// in : dx => spessore g/cm2 |
| 936 |
Float_t dEP=0.; |
Float_t dEP=0.; |
| 937 |
|
|
| 938 |
// gamma // Mass = A * MassP; /*in Mev/c2*/ |
// gamma // Mass = A * MassP; /*in Mev/c2*/ |
| 939 |
gam = (*Etot)/(elem[*Z-1]*MassP); // E = gamma M c2 |
gam = (*Etot)/(elem[*Z-1][*isotope]*MassP); // E = gamma M c2 |
| 940 |
|
|
| 941 |
|
|
| 942 |
Bet = sqrt((SQ(gam) -1.)/SQ(gam)); |
Bet = sqrt((SQ(gam) -1.)/SQ(gam)); |
| 962 |
|
|
| 963 |
|
|
| 964 |
|
|
| 965 |
void CaloBragg::Enetrack(Int_t* Z, Float_t* E0, Float_t* primo,Float_t* ultimo, Float_t out[]){ |
void CaloBragg::Enetrack(Int_t* Z, Int_t* isotope, Float_t* E0, Float_t* primo,Float_t* ultimo, Float_t out[]){ |
| 966 |
|
|
| 967 |
//calcola energia rilasciata sulla traccia (usa ELOSS) |
//calcola energia rilasciata sulla traccia (usa ELOSS) |
| 968 |
// in : Z =>carica |
// in : Z =>carica |
| 980 |
//azzero energia rilasciata sui piani |
//azzero energia rilasciata sui piani |
| 981 |
memset(out, 0, 2*NPLA*sizeof(Float_t)); |
memset(out, 0, 2*NPLA*sizeof(Float_t)); |
| 982 |
|
|
| 983 |
Float_t Massa = (elem[(*Z)-1] * MassP); |
Float_t Massa = (elem[(*Z)-1][*isotope] * MassP); |
| 984 |
|
|
| 985 |
for( Int_t ipla=((int)(*primo)); ipla<= ((int)(*ultimo)); ipla++){ |
for( Int_t ipla=((int)(*primo)); ipla<= ((int)(*ultimo)); ipla++){ |
| 986 |
dE=0.; |
dE=0.; |
| 987 |
//spessore silicio corretto x inclinazione, z, energia, out:rilascio |
//spessore silicio corretto x inclinazione, z, energia, out:rilascio |
| 988 |
ELOSS(&spessore[0], Z, &Ezero, &dE, ISi);//spessore in g/cm2!! |
ELOSS(&spessore[0], Z , isotope , &Ezero, &dE, ISi);//spessore in g/cm2!! |
| 989 |
if((Ezero-dE) <= Massa){//se l'energia depositata e' maggiore dell'energia della perticella stop |
|
| 990 |
|
if(dE!=dE) return; //controlla che non sia un NaN |
| 991 |
|
|
| 992 |
|
if((Ezero-dE) <= Massa){//se l'energia depositata e' maggiore dell'energia della perticella stop |
| 993 |
out[ipla] = Ezero - Massa; //MeV |
out[ipla] = Ezero - Massa; //MeV |
| 994 |
return; |
return; |
| 995 |
|
|
| 996 |
}else{ |
}else{ |
| 997 |
out[ipla] = dE; //MeV |
out[ipla] = dE; //MeV |
| 998 |
Ezero = Ezero - dE;//energia residua |
Ezero = Ezero - dE;//energia residua |
| 999 |
if ( debug ) printf(" zompa %i out %f dE %f ezero %f \n",ipla,out[ipla],dE,Ezero); |
// if ( debug ) printf(" zompa %i out %f dE %f ezero %f \n",ipla,out[ipla],dE,Ezero); |
| 1000 |
}; |
}; |
| 1001 |
//se sono su un piano Y (tutti i pari) dopo c'e' il tungsteno |
//se sono su un piano Y (tutti i pari) dopo c'e' il tungsteno |
| 1002 |
if(ipla%2 == 0){ |
if(ipla%2 == 0){ |
| 1011 |
sp = spessore[0]*spessore[1]; //((gcm2Si)*(WinSi))//spessore attraversato in g/cm2 |
sp = spessore[0]*spessore[1]; //((gcm2Si)*(WinSi))//spessore attraversato in g/cm2 |
| 1012 |
} |
} |
| 1013 |
// printf(" sp %f II %f \n",sp,II); |
// printf(" sp %f II %f \n",sp,II); |
| 1014 |
ELOSS(&sp, Z, &Ezero, &dE,II); |
ELOSS(&sp, Z, isotope , &Ezero, &dE,II); |
| 1015 |
if((Ezero-dE) <= Massa){//se l'energia depositata e' maggiore dell'energia della perticella stop |
if((Ezero-dE) <= Massa){//se l'energia depositata e' maggiore dell'energia della perticella stop |
| 1016 |
return; |
return; |
| 1017 |
}else{ |
}else{ |
| 1086 |
//meno peso ai piani con rilasci maggiori di 1000 MIP |
//meno peso ai piani con rilasci maggiori di 1000 MIP |
| 1087 |
// if(calorimetro[ipla][1] > 1000) wi=0.5; |
// if(calorimetro[ipla][1] > 1000) wi=0.5; |
| 1088 |
if(calorimetro[ipla][1] > 1200.) wi=0.5; |
if(calorimetro[ipla][1] > 1200.) wi=0.5; |
| 1089 |
|
if(debug) printf("chiquadro start \n "); |
| 1090 |
Float_t arg = w*wi*(dE[ipla] - (calorimetro[ipla][1] * MIP)); |
Float_t arg = w*wi*(dE[ipla] - (calorimetro[ipla][1] * MIP)); |
| 1091 |
|
|
| 1092 |
sum += SQ(arg); // w*wi*(dEpiani[p][v]-(eplane[p][v]*MIP))));//( dEpiani[p][v] - (eplane[p][v]*MIP)); |
sum += SQ(arg); // w*wi*(dEpiani[p][v]-(eplane[p][v]*MIP))));//( dEpiani[p][v] - (eplane[p][v]*MIP)); |
| 1127 |
//out: array[4]=> chi2,Zbest,Ebest,piani saltati nel chi2 |
//out: array[4]=> chi2,Zbest,Ebest,piani saltati nel chi2 |
| 1128 |
// |
// |
| 1129 |
|
|
| 1130 |
|
//printf("entrato"); |
| 1131 |
memset(dEplan,0,2*NPLA*sizeof(Float_t)); |
memset(dEplan,0,2*NPLA*sizeof(Float_t)); |
| 1132 |
|
|
| 1133 |
Int_t Z = 0;// z iniziale |
Int_t Z = 0;// z iniziale |
| 1134 |
|
|
| 1135 |
|
Int_t isotope=0; |
| 1136 |
|
|
| 1137 |
Float_t Massa = 0.; |
Float_t Massa = 0.; |
| 1138 |
|
|
| 1146 |
Int_t max=32;//max z di cui so la massa :P |
Int_t max=32;//max z di cui so la massa :P |
| 1147 |
if((Zlimite)<=31) max=(int)(Zlimite) + 1; |
if((Zlimite)<=31) max=(int)(Zlimite) + 1; |
| 1148 |
|
|
| 1149 |
|
if(debug) printf("loopze inizio max %d \n",max); |
| 1150 |
if ( fzeta > 0. ){ |
if ( fzeta > 0. ){ |
| 1151 |
zmin = fzeta; |
zmin = fzeta; |
| 1152 |
max = fzeta+1; |
max = fzeta+1; |
| 1154 |
|
|
| 1155 |
Int_t colmax=32; |
Int_t colmax=32; |
| 1156 |
Int_t rowmax=3000; |
Int_t rowmax=3000; |
| 1157 |
|
Int_t isomax=7; |
| 1158 |
|
|
| 1159 |
Float_t matrixchi2[colmax][rowmax][3]; |
Float_t matrixchi2[colmax][isomax][rowmax][3]; |
| 1160 |
memset(matrixchi2, 0, colmax*rowmax*3*sizeof(Float_t)); |
memset(matrixchi2, 0, colmax*isomax*rowmax*3*sizeof(Float_t)); |
| 1161 |
|
|
| 1162 |
Int_t imin = 1-nostep/2; |
Int_t imin = 1-nostep/2; |
| 1163 |
Int_t imax = nostep/2; |
Int_t imax = nostep/2; |
| 1166 |
for(Int_t inucl=zmin; inucl<max; inucl++){ |
for(Int_t inucl=zmin; inucl<max; inucl++){ |
| 1167 |
|
|
| 1168 |
Z= inucl; |
Z= inucl; |
| 1169 |
|
|
| 1170 |
Massa = elem[inucl-1]*MassP; |
//loop isotopi |
| 1171 |
|
while ( elem[inucl-1][isotope] > 0. ){ |
| 1172 |
|
|
| 1173 |
|
if( fiso != -1 ){ |
| 1174 |
|
isotope=fiso; |
| 1175 |
|
if(debug) printf("In Loopze - Isotope N %d",isotope); |
| 1176 |
|
} |
| 1177 |
|
Massa = elem[inucl-1][isotope]*MassP; |
| 1178 |
|
|
| 1179 |
//loop energia |
//loop energia |
| 1180 |
Int_t iene2 = 0; |
Int_t iene2 = -1; |
| 1181 |
|
|
| 1182 |
// for(Int_t iene= 0; iene<1000; iene++){// da non cambiare in base a Stepint altrimenti cambia la matrice bestchi2!!!cosi' non raggiungo mai integrale!!!!! mettere <=?? |
// for(Int_t iene= 0; iene<1000; iene++){// da non cambiare in base a Stepint altrimenti cambia la matrice bestchi2!!!cosi' non raggiungo mai integrale!!!!! mettere <=?? |
| 1183 |
|
|
| 1184 |
for(Int_t iene= imin; iene<imax; iene++){// da non cambiare in base a Stepint altrimenti cambia la matrice bestchi2!!!cosi' non raggiungo mai integrale!!!!! mettere <=?? |
for(Int_t iene= imin; iene<imax; iene++){// da non cambiare in base a Stepint altrimenti cambia la matrice bestchi2!!!cosi' non raggiungo mai integrale!!!!! mettere <=?? |
| 1185 |
|
iene2++; |
| 1186 |
|
energia= Massa + (E0)+ iene*Stepint;//gli do un'energia totale (momento) massa+energia cinetica, aumentando la cinetica.. |
| 1187 |
|
|
| 1188 |
|
|
| 1189 |
|
if( fene > 0. ) energia=fene; //forza l'energia |
| 1190 |
|
if (debug) printf("loopze energia %f, z %d, isotopo %d ,iene %d\n",energia,Z,isotope,iene); |
| 1191 |
|
// printf(" energia %f , forzata %f \n",energia,fene); |
| 1192 |
|
Enetrack(&Z, &isotope, &energia, &estremi[0][0],&estremi[1][0], dEplan);//calcola rilascio energetico sui piani |
| 1193 |
|
|
| 1194 |
iene2++; |
chiquadro(dEplan,chi2); //calcolo chi2 |
| 1195 |
energia= Massa + (E0)+ iene*Stepint;//gli do un'energia totale (momento) massa+energia cinetica, aumentando la cinetica.. |
if (debug) printf("loopze chi %f \n",chi2[0]); |
| 1196 |
|
if(debug && TMath::Finite(chi2[0])==1 && (TMath::IsNaN(chi2[0])!=1) ) printf("loopze fin mat %f \n",chi2[0]); |
|
Enetrack(&Z, &energia, &estremi[0][0],&estremi[1][0], dEplan);//calcola rilascio energetico sui piani |
|
|
//calcolo chi2 |
|
|
chiquadro(dEplan,chi2); |
|
|
|
|
| 1197 |
// printf(" last deplan from: Z = %i iene %i energia %f chi2 %f \n",inucl,iene,energia,chi2[0]); |
// printf(" last deplan from: Z = %i iene %i energia %f chi2 %f \n",inucl,iene,energia,chi2[0]); |
| 1198 |
|
if( (chi2[1] != 79.) ){//salto quelli che frammentano |
| 1199 |
if( (chi2[1] != 79.) ){//salto quelli che frammentano |
matrixchi2[inucl][isotope][iene2][0]=chi2[0];//valore chi2 per questo z a questa energia |
| 1200 |
matrixchi2[inucl][iene2][0]=chi2[0];//valore chi2 per questo z a questa energia |
matrixchi2[inucl][isotope][iene2][1]=energia;//energia per questo chi2 |
| 1201 |
matrixchi2[inucl][iene2][1]=energia;//energia per questo chi2 |
matrixchi2[inucl][isotope][iene2][2]=chi2[2];//piani saltati nel chi2 |
| 1202 |
matrixchi2[inucl][iene2][2]=chi2[2];//piani saltati nel chi2 |
if( fene > 0. ) break; |
| 1203 |
} else { |
} else { |
| 1204 |
matrixchi2[inucl][iene2][0]=1000.;//valore chi2 per questo z a questa energia |
matrixchi2[inucl][isotope][iene2][0]=numeric_limits<Float_t>::max();//valore chi2 per questo z a questa energia |
| 1205 |
matrixchi2[inucl][iene2][1]=1000.;//energia per questo chi2 |
matrixchi2[inucl][isotope][iene2][1]=numeric_limits<Float_t>::max();//energia per questo chi2 |
| 1206 |
matrixchi2[inucl][iene2][2]=1000.;//piani saltati nel chi2 |
matrixchi2[inucl][isotope][iene2][2]=numeric_limits<Float_t>::max();//piani saltati nel chi2 |
| 1207 |
break; |
break; |
| 1208 |
} |
} |
|
}//fine loop energia |
|
| 1209 |
|
|
| 1210 |
|
}//fine loop energia |
| 1211 |
|
|
| 1212 |
|
if( fiso != -1 ){ |
| 1213 |
|
if(debug) printf("exited form isotopes loop"); |
| 1214 |
|
break; |
| 1215 |
|
} |
| 1216 |
|
|
| 1217 |
|
isotope++; //incremento il contatore isotopi |
| 1218 |
|
}//fine loop isotopi |
| 1219 |
|
isotope=0; //riazzero il contatore isotopi |
| 1220 |
|
|
| 1221 |
};//fine loop z |
}//fine loop z |
| 1222 |
|
|
| 1223 |
|
isotope=0;//non dovrebbe servire |
| 1224 |
|
|
| 1225 |
//Emi |
//Emi |
| 1226 |
for (Int_t nu=zmin; nu<max; nu++){ |
for (Int_t nu=zmin; nu<max; nu++){ |
| 1227 |
for (Int_t en=0; en<nostep; en++){ |
if( fiso != -1 ){ |
| 1228 |
if((matrixchi2[nu][en][0]<bestchi2[0]) && (matrixchi2[nu][en][0] >0.)){ |
isotope=fiso; |
| 1229 |
bestchi2[0]= matrixchi2[nu][en][0];// chi2 |
if(debug) printf("In Loopze EMI - Isotope N %d",isotope); |
| 1230 |
bestchi2[1]= (Float_t)nu; // z |
} |
| 1231 |
bestchi2[2]= matrixchi2[nu][en][1];//energia; |
while(elem[nu-1][isotope]> 0.){ |
| 1232 |
bestchi2[3]= matrixchi2[nu][en][2];// totale piani saltati |
for (Int_t en=0; en<nostep; en++){ |
| 1233 |
} |
if((matrixchi2[nu][isotope][en][0]<bestchi2[0]) && (matrixchi2[nu][isotope][en][0] >0.)){ |
| 1234 |
} |
bestchi2[0]= matrixchi2[nu][isotope][en][0];// chi2 |
| 1235 |
|
bestchi2[1]= (Float_t)nu; // z |
| 1236 |
|
bestchi2[2]= matrixchi2[nu][isotope][en][1];//energia; |
| 1237 |
|
bestchi2[3]= matrixchi2[nu][isotope][en][2];// totale piani saltati |
| 1238 |
|
bestchi2[4]= (Float_t)isotope; //isotopo |
| 1239 |
|
} |
| 1240 |
|
} |
| 1241 |
|
|
| 1242 |
|
if( fiso != -1 ){ |
| 1243 |
|
if(debug) printf("exited form isotopes loop"); |
| 1244 |
|
break; |
| 1245 |
} |
} |
| 1246 |
|
|
| 1247 |
|
isotope++; |
| 1248 |
|
} |
| 1249 |
|
isotope=0; |
| 1250 |
|
} |
| 1251 |
|
|
| 1252 |
};//endloopze |
};//endloopze |
| 1253 |
|
|
| 1371 |
|
|
| 1372 |
Float_t zmin=1.; |
Float_t zmin=1.; |
| 1373 |
Float_t zmax=32.; |
Float_t zmax=32.; |
| 1374 |
Float_t bestchitemp[4] = {0,0,0,0}; |
Float_t bestchitemp[5] = {0,0,0,0,0}; |
| 1375 |
|
|
| 1376 |
bestchi2[0]=10000.; |
bestchi2[0]=numeric_limits<Float_t>::max(); |
| 1377 |
bestchi2[1]=0.; |
bestchi2[1]=0.; |
| 1378 |
bestchi2[2]=0.; |
bestchi2[2]=0.; |
| 1379 |
bestchi2[3]=0.; |
bestchi2[3]=0.; |
| 1380 |
Float_t zero=0.; |
bestchi2[4]=0.; |
| 1381 |
|
// Float_t zero=0.; |
| 1382 |
//------------primo loop ---------------------- |
//------------primo loop ---------------------- |
| 1383 |
// energia ezero, zstart zstop |
// energia ezero, zstart zstop |
| 1384 |
// loopze(Integrale,zero,zmin,zmax); |
// loopze(Integrale,zero,zmin,zmax); |
| 1387 |
loopze(Integrale*1.2/500.,Integrale/1000.,zmin,zmax,200); |
loopze(Integrale*1.2/500.,Integrale/1000.,zmin,zmax,200); |
| 1388 |
|
|
| 1389 |
// loopze(Integrale*2.,Integrale/100.,zmin,zmax); |
// loopze(Integrale*2.,Integrale/100.,zmin,zmax); |
| 1390 |
if ( debug ) printf(" Integrale %f , outene %f \n",Integrale,bestchi2[2]); |
if ( debug) printf("Zdaloop start Integrale %f , outene %f \n",Integrale,bestchi2[2]); |
| 1391 |
|
|
| 1392 |
//------------secondo loop ---------------------- |
//------------secondo loop ---------------------- |
| 1393 |
for(Int_t i=0;i<4;i++) bestchitemp[i]=bestchi2[i]; |
for(Int_t i=0;i<5;i++) bestchitemp[i]=bestchi2[i]; |
| 1394 |
bestchi2[0] = 10000.; |
bestchi2[0]=numeric_limits<Float_t>::max(); |
| 1395 |
bestchi2[1] = 0.; |
bestchi2[1] = 0.; |
| 1396 |
bestchi2[2] = 0.; |
bestchi2[2] = 0.; |
| 1397 |
bestchi2[3] = 0.;//riazzero |
bestchi2[3] = 0.; |
| 1398 |
|
bestchi2[4] = 0.;//riazzero |
| 1399 |
|
|
| 1400 |
Float_t step = bestchitemp[2];// |
Float_t step = bestchitemp[2];// |
| 1401 |
zero=0.; // qualsiasi altro valore peggiora le cose |
// zero=0.; // qualsiasi altro valore peggiora le cose |
| 1402 |
// zmin=zmax=bestchitemp[1]; |
// zmin=zmax=bestchitemp[1]; |
| 1403 |
zmin=bestchitemp[1]-1; |
zmin=bestchitemp[1]-1; |
| 1404 |
zmax=bestchitemp[1]+1; |
zmax=bestchitemp[1]+1; |
| 1405 |
// loopze(step,zero,zmin,zmax); // |
//loopze(step,zero,zmin,zmax); // |
| 1406 |
|
|
| 1407 |
//-> loopze(step,step/2.,zmin,zmax,200); // |
//-> loopze(step,step/2.,zmin,zmax,200); // |
| 1408 |
loopze(step,step/2.,zmin,zmax,500); // |
loopze(step,step/2.,zmin,zmax,500); // |
| 1409 |
|
|
| 1410 |
if ( debug ) printf(" Integrale2 %f , outene %f step %f \n",Integrale,bestchi2[2],step); |
//step = bestchitemp[2];// |
| 1411 |
|
|
| 1412 |
|
//loopze(step/2,3*step/4.,zmin,zmax,500); // |
| 1413 |
|
|
| 1414 |
|
if ( debug ) printf("Zdaloop Integrale2 %f , outene %f step %f \n",Integrale,bestchi2[2],step); |
| 1415 |
|
|
|
|
|
| 1416 |
//chi2,z,Etot,Pskip |
//chi2,z,Etot,Pskip |
| 1417 |
lpchi2=bestchi2[0]; |
lpchi2=bestchi2[0]; |
| 1418 |
lpz=bestchi2[1]; |
lpz=bestchi2[1]; |
| 1419 |
lpetot=bestchi2[2]; |
lpetot=bestchi2[2]; |
| 1420 |
lppskip=bestchi2[3]; |
lppskip=bestchi2[3]; |
| 1421 |
|
lpisotope=bestchi2[4]; |
| 1422 |
};//endZdaloop |
};//endZdaloop |
| 1423 |
|
|
| 1424 |
|
|