25 |
Float_t rig[18]; |
Float_t rig[18]; |
26 |
Float_t rmean[17]; |
Float_t rmean[17]; |
27 |
Int_t ntot[17]; |
Int_t ntot[17]; |
28 |
|
Int_t MDIM = 8213; |
29 |
|
//Int_t MDIM = 4128; |
30 |
//Float_t qqplane[17][43]; |
//Float_t qqplane[17][43]; |
31 |
//Int_t nnqplane[17][43]; |
//Int_t nnqplane[17][43]; |
32 |
|
|
35 |
TMatrixD *matrix[17]; |
TMatrixD *matrix[17]; |
36 |
TMatrixD *nmat[17]; |
TMatrixD *nmat[17]; |
37 |
|
|
38 |
TMatrixD *fqplane[17]; |
TMatrixD *fqplane; |
39 |
TMatrixD *fnqplane[17]; |
TMatrixD *fnqplane; |
40 |
TMatrixF *fmatrix[17]; |
//TMatrixD *fqplane[17]; |
41 |
|
//TMatrixD *fnqplane[17]; |
42 |
|
//TMatrixF *fmatrix[17]; |
43 |
|
//TMatrixF *fnmat[17]; |
44 |
|
//TMatrixD *fmatrix; |
45 |
|
//TMatrixD *fnmat; |
46 |
|
TMatrixF *fmatrix; |
47 |
|
//TMatrixF *fnmat; |
48 |
TMatrixF *fnmat[17]; |
TMatrixF *fnmat[17]; |
49 |
|
//Int_t finmat[43][191]; |
50 |
|
|
51 |
//=============================================================================== |
//=============================================================================== |
52 |
bool Select( PamLevel2* event ){ |
bool Select( PamLevel2* event ){ |
250 |
if ( view == 1 && (plane >0 || strip > ti0) ) break; |
if ( view == 1 && (plane >0 || strip > ti0) ) break; |
251 |
}; |
}; |
252 |
// if ( event->GetCaloLevel1()->qtotpl(0) > 7. ) return false; |
// if ( event->GetCaloLevel1()->qtotpl(0) > 7. ) return false; |
253 |
|
|
254 |
|
// if ( rigidity > 2.2 || rigidity < 1.5 ) return false; |
255 |
|
// printf(" rig %f CRIG %i SRIG %i \n",rigidity,CRIG,SRIG); |
256 |
// |
// |
257 |
return true; |
return true; |
258 |
} |
} |
303 |
// |
// |
304 |
memset(rmean, 0, 17*sizeof(Float_t)); |
memset(rmean, 0, 17*sizeof(Float_t)); |
305 |
memset(ntot, 0, 17*sizeof(Int_t)); |
memset(ntot, 0, 17*sizeof(Int_t)); |
306 |
|
// memset(finmat, 0, 43*191*sizeof(Int_t)); |
307 |
// |
// |
308 |
for (Int_t i=0; i < 17 ; i++){ |
for (Int_t i=0; i < 17 ; i++){ |
309 |
|
//for (Int_t i=3; i < 4 ; i++){ |
310 |
if ( !FULL ){ |
if ( !FULL ){ |
311 |
matrix[i] = new TMatrixD(43,43); |
matrix[i] = new TMatrixD(43,43); |
312 |
qplane[i] = new TArrayF(43); |
qplane[i] = new TArrayF(43); |
314 |
nmat[i] = new TMatrixD(43,43); |
nmat[i] = new TMatrixD(43,43); |
315 |
} else { |
} else { |
316 |
if ( MATRIX ){ |
if ( MATRIX ){ |
317 |
fmatrix[i] = new TMatrixF(4128,4128); |
// fmatrix = new TMatrixF(4128,4128); |
318 |
fnmat[i] = new TMatrixF(4128,4128); |
// fnmat = new TMatrixF(4128,4128); |
319 |
|
// fmatrix = new TMatrixF(8213,8213); |
320 |
|
// fnmat = new TMatrixF(8213,8213); |
321 |
|
fmatrix = new TMatrixF(MDIM,MDIM); |
322 |
|
// fnmat = new TMatrixF(MDIM,MDIM); |
323 |
|
fnmat[i] = new TMatrixF(43,191); |
324 |
|
// cf->WriteFullMatrix(fmatrix, i); |
325 |
|
// cf->WriteFullNMatrix(fnmat, i); |
326 |
|
// delete fmatrix; |
327 |
|
// delete fnmat; |
328 |
//fnmat[i] = new TMatrixI(8213,8213); |
//fnmat[i] = new TMatrixI(8213,8213); |
329 |
} else { |
} else { |
330 |
fqplane[i] = new TMatrixD(43,191); // 43 planes x 191 strip (= 1 + 95 x 2, one strip is the one transversed by the track that could be on the extreme right or left) |
fqplane = new TMatrixD(43,191); // 43 planes x 191 strip (= 1 + 95 x 2, one strip is the one transversed by the track that could be on the extreme right or left) |
331 |
fnqplane[i] = new TMatrixD(43,191);// |
fnqplane = new TMatrixD(43,191);// |
332 |
|
// |
333 |
|
cf->WriteFullMean(fqplane, i); |
334 |
|
cf->WriteFullNMean(fnqplane, i); |
335 |
|
delete fqplane; |
336 |
|
delete fnqplane; |
337 |
|
// |
338 |
}; |
}; |
339 |
}; |
}; |
340 |
}; |
}; |
394 |
// |
// |
395 |
// FULL CALORIMETER |
// FULL CALORIMETER |
396 |
// |
// |
397 |
|
fqplane = cf->LoadFullAverage(rbi); |
398 |
|
fnqplane = cf->LoadFullNAverage(rbi); |
399 |
CaloTrkVar *ct = L2->GetTrack(0)->GetCaloTrack(); |
CaloTrkVar *ct = L2->GetTrack(0)->GetCaloTrack(); |
400 |
// |
// |
401 |
Int_t nplane = 0; |
Int_t nplane = 0; |
419 |
// |
// |
420 |
for (Int_t k=0; k<191; k++){ |
for (Int_t k=0; k<191; k++){ |
421 |
mstrip = cd + k; |
mstrip = cd + k; |
422 |
if ( mstrip < (191-cs) ) (*fnqplane[rbi])[nplane][mstrip] += 1.; |
// if ( mstrip < (191-cs) ) (*fnqplane[rbi])[nplane][mstrip] += 1.; |
423 |
|
if ( mstrip < (191-cs) ) (*fnqplane)[nplane][mstrip] += 1.; |
424 |
}; |
}; |
425 |
}; |
}; |
426 |
}; |
}; |
439 |
// |
// |
440 |
mstrip = cd + strip; |
mstrip = cd + strip; |
441 |
// |
// |
442 |
(*fqplane[rbi])[nplane][mstrip] += mip; |
// (*fqplane[rbi])[nplane][mstrip] += mip; |
443 |
|
(*fqplane)[nplane][mstrip] += mip; |
444 |
// |
// |
445 |
}; |
}; |
446 |
|
// |
447 |
|
cf->WriteFullMean(fqplane, rbi); |
448 |
|
cf->WriteFullNMean(fnqplane, rbi); |
449 |
|
cf->UnLoadFullAverage(rbi); |
450 |
|
cf->UnLoadFullNAverage(rbi); |
451 |
|
delete fqplane; |
452 |
|
delete fnqplane; |
453 |
// |
// |
454 |
}; |
}; |
455 |
} |
} |
474 |
// |
// |
475 |
}; |
}; |
476 |
} else { |
} else { |
477 |
|
// |
478 |
for (Int_t i=0; i<nbin-1; i++){ |
for (Int_t i=0; i<nbin-1; i++){ |
479 |
|
fqplane = cf->LoadFullAverage(i); |
480 |
|
fnqplane = cf->LoadFullNAverage(i); |
481 |
if ( ntot[i] > 0 ) rmean[i] /= (Float_t)(ntot[i]); |
if ( ntot[i] > 0 ) rmean[i] /= (Float_t)(ntot[i]); |
482 |
// |
// |
483 |
for (Int_t j=0; j<43 ; j++){ |
for (Int_t j=0; j<43 ; j++){ |
484 |
for (Int_t k=0; k<191; k++){ |
for (Int_t k=0; k<191; k++){ |
485 |
if ( (*fnqplane[i])[j][k] > 0 ){ |
// if ( (*fnqplane[i])[j][k] > 0 ){ |
486 |
(*fqplane[i])[j][k] /= (Float_t)(*fnqplane[i])[j][k]; |
// (*fqplane[i])[j][k] /= (Float_t)(*fnqplane[i])[j][k]; |
487 |
|
// } else { |
488 |
|
// (*fqplane[i])[j][k] = 0.; |
489 |
|
// }; |
490 |
|
// printf(" BIN %i plane %i strip %i average energy %f qplane %f nqplane %f \n",i,j,k,rmean[i],(*fqplane[i])[j][k],(*fnqplane[i])[j][k]); |
491 |
|
if ( (*fnqplane)[j][k] > 0 ){ |
492 |
|
if ( (*fqplane)[j][k] == 0. ) (*fqplane)[j][k] = 0.7; |
493 |
|
(*fqplane)[j][k] /= (Float_t)(*fnqplane)[j][k]; |
494 |
} else { |
} else { |
495 |
(*fqplane[i])[j][k] = 0.; |
(*fqplane)[j][k] = 0.; |
496 |
}; |
}; |
497 |
printf(" BIN %i plane %i strip %i average energy %f qplane %f nqplane %f \n",i,j,k,rmean[i],(*fqplane[i])[j][k],(*fnqplane[i])[j][k]); |
printf(" BIN %i plane %i strip %i average energy %f qplane %f nqplane %f \n",i,j,k,rmean[i],(*fqplane)[j][k],(*fnqplane)[j][k]); |
498 |
}; |
}; |
499 |
}; |
}; |
500 |
|
cf->WriteFullMean(fqplane, i); |
501 |
|
cf->WriteFullNMean(fnqplane, i); |
502 |
|
cf->UnLoadFullAverage(i); |
503 |
|
cf->UnLoadFullNAverage(i); |
504 |
|
delete fqplane; |
505 |
|
delete fnqplane; |
506 |
}; |
}; |
507 |
for (Int_t i=0; i<nbin-1; i++){ |
// |
508 |
// |
// for (Int_t i=0; i<nbin-1; i++){ |
509 |
cf->WriteFullMean(fqplane[i], i); |
// // |
510 |
// |
// cf->WriteFullMean(fqplane[i], i); |
511 |
}; |
// // |
512 |
|
// }; |
513 |
}; |
}; |
514 |
// |
// |
515 |
cf->WriteNumBin(nbin); |
cf->WriteNumBin(nbin); |
540 |
}; |
}; |
541 |
}; |
}; |
542 |
// |
// |
543 |
|
if ( rbi != 3 ) return; |
544 |
if ( erig < rig[0] ) return; |
if ( erig < rig[0] ) return; |
545 |
if ( erig >= rig[nbin-1] ) return; |
if ( erig >= rig[nbin-1] ) return; |
546 |
// |
// |
583 |
// |
// |
584 |
// FULL CALORIMETER |
// FULL CALORIMETER |
585 |
// |
// |
586 |
|
printf(" Retrieve matrix %i IEV %i \n",rbi,iev); |
587 |
|
// fmatrix = cf->LoadFullMatrix(rbi); |
588 |
|
// cf->LoadFullMatrix(rbi,fmatrix); |
589 |
|
// fnmat = cf->LoadFullNMatrix(rbi); |
590 |
|
printf(" done \n"); |
591 |
|
printf(" start loop \n"); |
592 |
|
// |
593 |
CaloTrkVar *ct = L2->GetTrack(0)->GetCaloTrack(); |
CaloTrkVar *ct = L2->GetTrack(0)->GetCaloTrack(); |
594 |
// |
// |
595 |
Int_t nplane = 0; |
Int_t nplane = 0; |
598 |
Int_t strip = 0; |
Int_t strip = 0; |
599 |
Float_t mip = 0.; |
Float_t mip = 0.; |
600 |
// |
// |
601 |
Int_t mindgf = 48; |
// Int_t mindgf = 48; |
602 |
Int_t dgf = 143; |
// Int_t dgf = 143; |
603 |
|
Int_t mindgf = 0; |
604 |
|
Int_t dgf = 191; |
605 |
|
// Int_t mindgf = 85; |
606 |
|
// Int_t dgf = 115; |
607 |
Int_t cs = 0; |
Int_t cs = 0; |
608 |
Int_t cd = 0; |
Int_t cd = 0; |
609 |
Int_t mstrip = 0; |
Int_t mstrip = 0; |
628 |
// |
// |
629 |
}; |
}; |
630 |
// |
// |
631 |
Float_t mip1; |
Float_t mip1 = 1.; |
632 |
Int_t cs1; |
Int_t cs1; |
633 |
Int_t cd1; |
Int_t cd1; |
634 |
Float_t mip2; |
Float_t mip2 = 1.; |
635 |
Int_t cs2; |
Int_t cs2; |
636 |
Int_t cd2; |
Int_t cd2; |
637 |
Int_t mi = -1; |
Int_t mi = -1; |
647 |
Int_t mstrip2min = 0; |
Int_t mstrip2min = 0; |
648 |
Int_t mstrip2max = 0; |
Int_t mstrip2max = 0; |
649 |
// |
// |
650 |
|
Int_t toto = 0; |
651 |
|
// |
652 |
for (Int_t nplane1=0; nplane1<43; nplane1++){ |
for (Int_t nplane1=0; nplane1<43; nplane1++){ |
653 |
for (Int_t mstrip1=0; mstrip1<191; mstrip1++){ |
if ( nplane1 >= 37 ) nn1 = nplane1 + 1; |
654 |
|
vi1 = 1; |
655 |
|
if ( nn1%2 ) vi1 = 0; |
656 |
|
pl1 = (nn1 - 1 + vi1)/2; |
657 |
|
// |
658 |
|
cs1 = ct->tibar[pl1][vi1] - 1; // convertire nplane in pl1 e vi1 |
659 |
|
// |
660 |
|
cd1 = 95 - cs1; |
661 |
|
// |
662 |
|
mstrip1min = TMath::Max(mindgf,(cd1+0)); |
663 |
|
mstrip1max = TMath::Min(dgf,(cd1+95)) + 1; |
664 |
|
// |
665 |
|
if ( nplane1 == 0 || nplane1 == 42 ) printf(" pl %i mstrip1min %i mstrip1max %i \n",nplane1,mstrip1min,mstrip1max); |
666 |
|
// |
667 |
|
for (Int_t mstrip1=mstrip1min; mstrip1<mstrip1max; mstrip1++){ |
668 |
|
// printf(".\n"); |
669 |
|
// |
670 |
mj = -1; |
mj = -1; |
671 |
// |
// |
672 |
if ( nplane1 >= 37 ) nn1 = nplane1 + 1; |
mip1 = mipv[nplane1][mstrip1] - cf->GetFullAverageAt(nplane1,mstrip1,erig,rbi); |
673 |
vi1 = 1; |
// |
674 |
if ( nn1%2 ) vi1 = 0; |
mi = (nplane1 * 191) + mstrip1; |
675 |
pl1 = (nn1 - 1 + vi1)/2; |
// |
676 |
// |
// if ( mstrip1 > mstrip1min ) break; |
677 |
cs1 = ct->tibar[pl1][vi1] - 1; // convertire nplane in pl1 e vi1 |
// if ( mstrip1 > dgf ) break; |
678 |
// |
// if ( mstrip1 >= mindgf && mstrip1 <= dgf && mstrip1 >= mstrip1min && mstrip1 <= mstrip1max ){ |
679 |
cd1 = 95 - cs1; |
// |
680 |
// |
// finmat[nplane1][mstrip1]++; |
681 |
mstrip1min = cd1 + 0; |
(*fnmat[rbi])[nplane1][mstrip1] += 1.; |
682 |
mstrip1max = cd1 + 95; |
// |
683 |
mip1 = mipv[nplane1][mstrip1]; |
if ( mip1 != 0. ){ |
684 |
// |
// |
|
if ( mstrip1 >= mindgf && mstrip1 <= dgf ){ |
|
|
mi++; |
|
|
// |
|
685 |
for (Int_t nplane2=0; nplane2<43; nplane2++){ |
for (Int_t nplane2=0; nplane2<43; nplane2++){ |
686 |
for (Int_t mstrip2=0; mstrip2<191; mstrip2++){ |
// |
687 |
// |
if ( nplane2 >= 37 ) nn2 = nplane2 + 1; |
688 |
if ( nplane2 >= 37 ) nn2 = nplane2 + 1; |
vi2 = 1; |
689 |
vi2 = 1; |
if ( nn2%2 ) vi2 = 0; |
690 |
if ( nn2%2 ) vi2 = 0; |
pl1 = (nn2 - 1 + vi2)/2; |
691 |
pl1 = (nn2 - 1 + vi2)/2; |
// |
692 |
|
cs2 = ct->tibar[pl2][vi2] - 1; |
693 |
|
// |
694 |
|
cd2 = 95 - cs2; |
695 |
|
// |
696 |
|
// mstrip2min = cd2 + 0; |
697 |
|
// mstrip2max = cd2 + 95; |
698 |
|
mstrip2min = TMath::Max(mindgf,(cd2+0)); |
699 |
|
mstrip2max = TMath::Min(dgf,(cd2+95)) + 1; |
700 |
|
// |
701 |
|
if ( nplane1 == 0 && nplane2 == 0 && mstrip1==mstrip1min ) printf(" mstrip2min %i mstrip2max %i \n",mstrip2min,mstrip2max); |
702 |
|
// |
703 |
|
for (Int_t mstrip2=mstrip2min; mstrip2<mstrip2max; mstrip2++){ |
704 |
// |
// |
705 |
cs2 = ct->tibar[pl2][vi2] - 1; |
mip2 = mipv[nplane2][mstrip2] - cf->GetFullAverageAt(nplane2,mstrip2,erig,rbi); |
706 |
// |
// |
707 |
cd2 = 95 - cs2; |
if ( mip2 != 0. ){ |
708 |
// |
// |
709 |
mstrip2min = cd2 + 0; |
mj = (nplane2 * 191) + mstrip2; |
710 |
mstrip2max = cd2 + 95; |
// mj++; |
711 |
// |
// |
712 |
mip2 = mipv[nplane2][mstrip2]; |
// if ( mstrip2 > mstrip2min ) break; |
713 |
// |
// if ( mstrip2 > dgf ) break; |
714 |
if ( mstrip2 >= mindgf && mstrip2 <= dgf ){ |
// if ( mstrip2 >= mindgf && mstrip2 <= dgf && mstrip2 >= mstrip2min && mstrip2 <= mstrip2max ){ |
715 |
mj++; |
// if ( mstrip1 >= mstrip1min && mstrip1 <= mstrip1max && mstrip2 >= mstrip2min && mstrip2 <= mstrip2max){ |
716 |
(*fmatrix[rbi])[mi][mj] += (mipv[nplane1][mstrip1] - cf->GetFullAverageAt(nplane1,mstrip1,erig)) * (mipv[nplane2][mstrip2] - cf->GetFullAverageAt(nplane2,mstrip2,erig)); |
// (*fmatrix[rbi])[mi][mj] += (mipv[nplane1][mstrip1] - cf->GetFullAverageAt(nplane1,mstrip1,erig)) * (mipv[nplane2][mstrip2] - cf->GetFullAverageAt(nplane2,mstrip2,erig)); |
717 |
if ( mstrip1 >= mstrip1min && mstrip1 <= mstrip1max && mstrip2 >= mstrip2min && mstrip2 <= mstrip2max){ |
// (*fnmat[rbi])[mi][mj] += 1.; |
718 |
(*fnmat[rbi])[mi][mj] += 1.; |
(*fmatrix)[mi][mj] += (mip1 * mip2); // giusto |
719 |
}; |
toto++; |
720 |
|
// (*fmatrix)[mi][mj] += 1.; |
721 |
|
// cf->GetFullAverageAt(nplane1,mstrip1,erig,rbi); |
722 |
|
// cf->GetFullAverageAt(nplane2,mstrip2,erig,rbi); |
723 |
|
// (*fnmat)[mi][mj] += 1.; |
724 |
|
// }; |
725 |
}; |
}; |
726 |
}; |
}; |
727 |
}; |
}; |
729 |
}; |
}; |
730 |
}; |
}; |
731 |
// |
// |
732 |
|
printf(" toto = %i \n",toto); |
733 |
|
printf("\n done \n"); |
734 |
|
printf(" write matrix \n"); |
735 |
|
// cf->WriteFullMatrix(fmatrix, rbi); |
736 |
|
// cf->WriteFullNMatrix(fnmat, rbi); |
737 |
|
printf(" done \n"); |
738 |
|
printf(" unload matrix \n"); |
739 |
|
// cf->UnLoadFullMatrix(rbi); |
740 |
|
// cf->UnLoadFullNMatrix(rbi); |
741 |
|
printf(" done \n"); |
742 |
|
printf(" delete matrix \n"); |
743 |
|
// delete fmatrix; |
744 |
|
// delete fnmat; |
745 |
|
printf(" done \n"); |
746 |
}; |
}; |
747 |
} |
} |
748 |
|
|
798 |
// |
// |
799 |
// FULL |
// FULL |
800 |
// |
// |
801 |
for (Int_t i=0; i<nbin-1; i++){ |
// for (Int_t i=0; i<nbin-1; i++){ |
802 |
|
// for (Int_t i=3; i<5; i++){ |
803 |
|
for (Int_t i=3; i<4; i++){ |
804 |
// |
// |
805 |
// determine the average matrix |
// determine the average matrix |
806 |
// |
// |
807 |
for (Int_t ii=0; ii<4171; ii++){ |
// fmatrix = cf->LoadFullMatrix(i); |
808 |
for (Int_t j=0; j<4171; j++){ |
// fnmat = cf->LoadFullNMatrix(i); |
809 |
if ( (*fnmat[i])[ii][j] > 0. ){ |
// |
810 |
(*fmatrix[i])[ii][j] /= (*fnmat[i])[ii][j]; |
// for (Int_t ii=0; ii<MDIM; ii++){ |
811 |
} else { |
// for (Int_t j=0; j<MDIM; j++){ |
812 |
(*fmatrix[i])[ii][j] = 0.; |
// // if ( (*fnmat[i])[ii][j] > 0. ){ |
813 |
|
// // (*fmatrix[i])[ii][j] /= (*fnmat[i])[ii][j]; |
814 |
|
// // } else { |
815 |
|
// // (*fmatrix[i])[ii][j] = 0.; |
816 |
|
// // }; |
817 |
|
// if ( (*fnmat)[ii][j] > 0. ){ |
818 |
|
// (*fmatrix)[ii][j] /= (*fnmat)[ii][j]; |
819 |
|
// } else { |
820 |
|
// (*fmatrix)[ii][j] = 0.; |
821 |
|
// }; |
822 |
|
// }; |
823 |
|
// }; |
824 |
|
// |
825 |
|
Int_t i1 = -1; |
826 |
|
Int_t j1 = -1; |
827 |
|
Int_t nonzero = 0; |
828 |
|
Int_t nonzero1 = 0; |
829 |
|
for (Int_t ii=0; ii<43; ii++){ |
830 |
|
for (Int_t j=0; j<191; j++){ |
831 |
|
// if ( (*fnmat[i])[ii][j] > 0. ){ |
832 |
|
// (*fmatrix[i])[ii][j] /= (*fnmat[i])[ii][j]; |
833 |
|
// } else { |
834 |
|
// (*fmatrix[i])[ii][j] = 0.; |
835 |
|
// }; |
836 |
|
i1 = (ii * 191) + j; |
837 |
|
// j1 = -1; |
838 |
|
for (Int_t iij=0; iij<43; iij++){ |
839 |
|
for (Int_t jj=0; jj<191; jj++){ |
840 |
|
// |
841 |
|
j1 = (iij * 191) + jj; |
842 |
|
// j1++; |
843 |
|
// if ( finmat[ii][j] > 0 ){ |
844 |
|
// (*fmatrix)[i1][j1] /= finmat[ii][j]; |
845 |
|
if ( (*fnmat[i])[ii][j] == 0. || (*fmatrix)[i1][j1] == 0. || !((*fmatrix)[i1][j1] == (*fmatrix)[i1][j1]) ){ |
846 |
|
(*fmatrix)[i1][j1] = 0.; |
847 |
|
} else { |
848 |
|
(*fmatrix)[i1][j1] /= (*fnmat[i])[ii][j]; |
849 |
|
nonzero++; |
850 |
|
if ( i1 == 0 ) nonzero1++; |
851 |
|
}; |
852 |
|
}; |
853 |
}; |
}; |
854 |
}; |
}; |
855 |
}; |
}; |
856 |
// |
// |
857 |
cf->WriteFullMatrix(fmatrix[i],i); |
printf(" Matrix has %i non-zero elements \n",nonzero); |
858 |
|
printf(" Matrix has %i non-zero elements on the first row\n",nonzero1); |
859 |
|
// |
860 |
|
// Bool_t BAD = false; |
861 |
|
// for (Int_t ii=0; ii<43; ii++){ |
862 |
|
// for (Int_t j=0; j<191; j++){ |
863 |
|
// // |
864 |
|
// i1 = (ii * 191) + j; |
865 |
|
// // |
866 |
|
// for (Int_t iij=0; iij<43; iij++){ |
867 |
|
// for (Int_t jj=0; jj<191; jj++){ |
868 |
|
// // |
869 |
|
// j1 = (iij * 191) + jj; |
870 |
|
// // |
871 |
|
// // printf(" ROW %i COLUMN %i VALUE %f \n",i1,j1,(*fmatrix)[i1][j1]); |
872 |
|
// if ( (*fmatrix)[i1][j1] == 0. || !((*fmatrix)[i1][j1]==(*fmatrix)[i1][j1]) ){ |
873 |
|
// printf(" ROW %i COLUMN %i VALUE %f \n",i1,j1,(*fmatrix)[i1][j1]); |
874 |
|
// printf(" che schifo! \n"); |
875 |
|
// BAD = true; |
876 |
|
// }; |
877 |
|
// // |
878 |
|
// }; |
879 |
|
// }; |
880 |
|
// }; |
881 |
|
// }; |
882 |
|
// // |
883 |
|
// if ( BAD ) printf(" questa matrice fa cagare \n"); |
884 |
|
// |
885 |
|
// |
886 |
|
// cf->WriteFullMatrix(fmatrix[i],i); |
887 |
|
cf->WriteFullMatrix(fmatrix, i); |
888 |
|
// cf->WriteFullNMatrix(fnmat, i); |
889 |
|
cf->WriteFullNMatrix(fnmat[i], i); |
890 |
// |
// |
891 |
if ( fmatrix[i]->Determinant() == 0. ){ |
// if ( fmatrix[i]->Determinant() == 0. ){ |
892 |
|
if ( fmatrix->Determinant() == 0. ){ |
893 |
printf(" ERROR: the matrix at bin %i is singular, determinant = 0., it cannot be inverted! \n",i); |
printf(" ERROR: the matrix at bin %i is singular, determinant = 0., it cannot be inverted! \n",i); |
894 |
} else { |
} else { |
895 |
Double_t det = 0.; |
Double_t det = 0.; |
896 |
TMatrixF invmatrix = (TMatrixF)(fmatrix[i]->Invert(&det)); |
// TMatrixF invmatrix = (TMatrixF)(fmatrix[i]->Invert(&det)); |
897 |
|
// printf(" Bin %i determinant is %f \n",i,det); |
898 |
|
// cf->WriteInvertedFullMatrix((TMatrixF)invmatrix,i); |
899 |
|
TMatrixF invmatrix = (TMatrixF)(fmatrix->Invert(&det)); |
900 |
printf(" Bin %i determinant is %f \n",i,det); |
printf(" Bin %i determinant is %f \n",i,det); |
901 |
cf->WriteInvertedFullMatrix((TMatrixF)invmatrix,i); |
cf->WriteInvertedFullMatrix((TMatrixF)invmatrix,i); |
902 |
}; |
}; |
903 |
|
// |
904 |
|
cf->UnLoadFullMatrix(i); |
905 |
|
// cf->UnLoadFullNMatrix(i); |
906 |
|
delete fmatrix; |
907 |
|
// delete fnmat; |
908 |
|
// |
909 |
}; |
}; |
910 |
}; |
}; |
911 |
// |
// |