--- calo/flight/CaloFranzini/src/Calib.cpp 2008/01/21 10:24:10 1.7 +++ calo/flight/CaloFranzini/src/Calib.cpp 2008/01/25 15:09:07 1.8 @@ -58,6 +58,7 @@ PamTrack *track = event->GetTrack(0); if(!track)return false; + if ( SIMU ) return true; //------------------------------------------------------------------ // tracker pre-selection //------------------------------------------------------------------ @@ -264,16 +265,22 @@ // // //=============================================================== -void CreateHistos( PamLevel2* event , TString file){ +void CreateHistos( PamLevel2* event , TString file){ cf = new CaloFranzini(event); // + if ( FULL ){ + cf->CalculateLongTZeta(false); + cf->CalculateFullTZeta(); + }; + // if ( MATRIX ){ cf->UpdateMatrixFile(file.Data()); - cf->LoadBin(); if ( !FULL ){ + cf->LoadBin(); cf->LoadLong(); } else { + cf->LoadBin(true); cf->LoadFull(); }; } else { @@ -323,9 +330,11 @@ // fnmat = new TMatrixF(MDIM,MDIM); // fmatrix[i] = new TMatrixF(1849,1849); // fnmat[i] = new TMatrixF(43,43); - fmatrix[i] = new TMatrixD(1333,1333); + // fmatrix[i] = new TMatrixD(1333,1333); + // fnmat[i] = new TMatrixF(43,31); + fmatrix[i] = new TMatrixD(473,473); + fnmat[i] = new TMatrixF(43,11); // fmatrix[i]->SetTol(tol); - fnmat[i] = new TMatrixF(43,31); // cf->WriteFullMatrix(fmatrix, i); // cf->WriteFullNMatrix(fnmat, i); // delete fmatrix; @@ -337,8 +346,10 @@ // fqplane[i] = new TMatrixD(43,43); // 43 planes x 43 "strip", where 43 = 50 + 14 + 5 + 2 + 2 + 2 + 2 + 2 + 2 + 2 + 2 + 1 + 1 + 1 + 1 + 1 + 1 + 1 + 1 + 1 + 1 + [1] + ... // fnqplane[i] = new TMatrixD(43,43);// 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 // - fqplane[i] = new TMatrixD(43,31); // 43 planes x 43 "strip", where 43 = 50 + 14 + 6 + 5 + 3 + 3 + 3 + 2 + 2 + 2 + 1 + 1 + 1 + 1 + 1 + [1] + 1 + 1 + 1 + 1 + ... - fnqplane[i] = new TMatrixD(43,31);// 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 ... + // fqplane[i] = new TMatrixD(43,31); // 43 planes x 43 "strip", where 43 = 50 + 14 + 6 + 5 + 3 + 3 + 3 + 2 + 2 + 2 + 1 + 1 + 1 + 1 + 1 + [1] + 1 + 1 + 1 + 1 + ... + // fnqplane[i] = new TMatrixD(43,31);// 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 ... + fqplane[i] = new TMatrixD(43,11); // 43 planes x 11 "strip", where 43 = 84 + 6 + 3 + 1 + 1 +[1]+ 1 + 1 + 3 + 6 + 84 + fnqplane[i] = new TMatrixD(43,11);// 0 1 2 3 4 5 6 7 8 9 10 // // cf->WriteFullMean(fqplane, i); // cf->WriteFullNMean(fnqplane, i); @@ -501,7 +512,8 @@ for (Int_t j=0; j<43 ; j++){ // for (Int_t k=0; k<191; k++){ // for (Int_t k=0; k<43; k++){ - for (Int_t k=0; k<31; k++){ + // for (Int_t k=0; k<31; k++){ + for (Int_t k=0; k<11; k++){ // if ( (*fnqplane[i])[j][k] > 0 ){ // (*fqplane[i])[j][k] /= (Float_t)(*fnqplane[i])[j][k]; // } else { @@ -538,9 +550,15 @@ cf->WriteRigBin(rigbin); // TArrayF *rmeanbin = new TArrayF(17, rmean); - TFile *file = cf->GetFile(); - file->cd(); - file->WriteObject(&(*rmeanbin), "binrigmean"); + if ( FULL ){ + TFile *file = cf->GetFullFile(); + file->cd(); + file->WriteObject(&(*rmeanbin), "binrigmean"); + } else { + TFile *file = cf->GetFile(); + file->cd(); + file->WriteObject(&(*rmeanbin), "binrigmean"); + }; // // } @@ -627,15 +645,18 @@ // Int_t mindgf = 84; // Int_t dgf = 106; Int_t mindgf = 0; - Int_t dgf = 43; + // Int_t dgf = 43; + Int_t dgf = 11; Int_t cs = 0; Int_t cd = 0; Int_t mstrip = 0; // // Float_t mipv[43][43]; // memset(mipv,0,43*43*sizeof(Float_t)); - Float_t mipv[43][31]; - memset(mipv,0,43*31*sizeof(Float_t)); + // Float_t mipv[43][31]; + // memset(mipv,0,43*31*sizeof(Float_t)); + Float_t mipv[43][11]; + memset(mipv,0,43*11*sizeof(Float_t)); // for ( Int_t i=0; iGetCaloLevel1()->istrip; i++ ){ // @@ -676,6 +697,17 @@ // Int_t toto = 0; // + // + // + Int_t therigb = rbi; + // + if ( erig < rig[0] ){ + therigb = 0; + }; + if ( erig >= rig[nbin-2] ){ + therigb = nbin - 3; + }; + // for (Int_t nplane1=0; nplane1<43; nplane1++){ if ( nplane1 >= 37 ) nn1 = nplane1 + 1; vi1 = 1; @@ -700,11 +732,12 @@ // mj = -1; // - mip1 = mipv[nplane1][mstrip1] - cf->GetFullAverageAt(nplane1,mstrip1,erig,rbi); + mip1 = mipv[nplane1][mstrip1] - cf->GetFullAverageAt(nplane1,mstrip1,erig,therigb); // // mi = (nplane1 * 191) + mstrip1; // mi = (nplane1 * 43) + mstrip1; - mi = (nplane1 * 31) + mstrip1; + // mi = (nplane1 * 31) + mstrip1; + mi = (nplane1 * 11) + mstrip1; // // if ( mstrip1 > mstrip1min ) break; // if ( mstrip1 > dgf ) break; @@ -737,20 +770,35 @@ // for (Int_t mstrip2=mstrip2min; mstrip2GetFullAverageAt(nplane2,mstrip2,erig,rbi); + mip2 = mipv[nplane2][mstrip2] - cf->GetFullAverageAt(nplane2,mstrip2,erig,therigb); // if ( mip2 != 0. ){ // // mj = (nplane2 * 191) + mstrip2; // mj = (nplane2 * 43) + mstrip2; // mj = (nplane2 * 31) + mstrip2; - Int_t sh = -15 + nplane1; - if ( sh > 15 ) sh -= 31*nplane1; - // - mj = (nplane2 * 31) + mstrip2 + sh; +// Int_t sh = -15 + nplane1; +// if ( sh > 15 ) sh -= 31*nplane1; +// // +// mj = (nplane2 * 31) + mstrip2 + sh; +// // +// if ( mj < 0 ) mj += 1333; +// if ( mj >= 1333 ) mj -= 1333; +// +// +// +// +// Int_t sh = -5 + nplane1; +// if ( sh > 5 ) sh -= 11*nplane1; +// // +// mj = (nplane2 * 11) + mstrip2 + sh; +// // +// if ( mj < 0 ) mj += 473; +// if ( mj >= 473 ) mj -= 473; +// // + // + mj = (nplane2 * 11) + mstrip2; // - if ( mj < 0 ) mj += 1333; - if ( mj >= 1333 ) mj -= 1333; // printf(" mi %i mj %i sh %i \n",mi,mj,sh); // // mj++; @@ -885,7 +933,8 @@ for (Int_t ii=0; ii<43; ii++){ // for (Int_t j=0; j<191; j++){ // for (Int_t j=0; j<43; j++){ - for (Int_t j=0; j<31; j++){ + // for (Int_t j=0; j<31; j++){ + for (Int_t j=0; j<11; j++){ // if ( (*fnmat[i])[ii][j] > 0. ){ // (*fmatrix[i])[ii][j] /= (*fnmat[i])[ii][j]; // } else { @@ -893,24 +942,36 @@ // }; // i1 = (ii * 191) + j; // i1 = (ii * 43) + j; - i1 = (ii * 31) + j; + // i1 = (ii * 31) + j; + i1 = (ii * 11) + j; // j1 = -1; for (Int_t iij=0; iij<43; iij++){ // for (Int_t jj=0; jj<191; jj++){ // for (Int_t jj=0; jj<43; jj++){ - for (Int_t jj=0; jj<31; jj++){ + // for (Int_t jj=0; jj<31; jj++){ + for (Int_t jj=0; jj<11; jj++){ // // j1 = (iij * 191) + jj; // j1 = (iij * 43) + jj; - Int_t sh = -15 + ii; - if ( sh > 15 ) sh -= 31*ii; +// Int_t sh = -15 + ii; +// if ( sh > 15 ) sh -= 31*ii; +// // +// j1 = (iij * 31) + jj + sh; +// // +// if ( j1 < 0 ) j1 += 1333; +// if ( j1 >= 1333 ) j1 -= 1333; // - j1 = (iij * 31) + jj + sh; // - if ( j1 < 0 ) j1 += 1333; - if ( j1 >= 1333 ) j1 -= 1333; - -// j1 = (iij * 31) + jj; + // +// Int_t sh = -5 + ii; +// if ( sh > 5 ) sh -= 11*ii; +// // +// j1 = (iij * 11) + jj + sh; +// // +// if ( j1 < 0 ) j1 += 473; +// if ( j1 >= 473 ) j1 -= 473; +// + j1 = (iij * 11) + jj; // j1++; // if ( finmat[ii][j] > 0 ){ // (*fmatrix)[i1][j1] /= finmat[ii][j]; @@ -1001,10 +1062,10 @@ // if ( BAD ) printf(" questa matrice fa cagare \n"); // // - cf->WriteFullMatrix(fmatrix[i],i); // cf->WriteFullMatrix(fmatrix, i); // cf->WriteFullNMatrix(fnmat, i); - cf->WriteFullNMatrix(fnmat[i], i); + //cf->WriteFullMatrix(fmatrix[i],i);// OK + // cf->WriteFullNMatrix(fnmat[i], i); //OK // // TDecompSVD svd(*fmatrix[i]); // Bool_t ok = svd.Decompose(); @@ -1014,6 +1075,8 @@ if ( fmatrix[i]->Determinant() == zero ){ //if ( fmatrix->Determinant() == 0. ){ printf(" ERROR: the matrix at bin %i is singular, determinant = 0., it cannot be inverted! \n",i); + cf->WriteFullMatrix(fmatrix[i],i);// OK + cf->WriteFullNMatrix(fnmat[i], i); //OK } else { // }; // if ( i == 3 ){