00001 #include <iostream>
00002 #include <fstream>
00003 #include <vector>
00004 #include <string>
00005 #include <cmath>
00006
00007 #include <FWCore/Framework/interface/Frameworkfwd.h>
00008 #include <FWCore/Framework/interface/MakerMacros.h>
00009 #include "FWCore/MessageLogger/interface/MessageLogger.h"
00010 #include "FWCore/Framework/interface/EDAnalyzer.h"
00011 #include "FWCore/Framework/interface/Event.h"
00012 #include "DataFormats/Common/interface/Handle.h"
00013 #include "FWCore/ParameterSet/interface/ParameterSet.h"
00014 #include "FWCore/Framework/interface/EventSetup.h"
00015 #include "FWCore/Framework/interface/ESHandle.h"
00016 #include "DataFormats/CSCDigi/interface/CSCStripDigi.h"
00017 #include "DataFormats/CSCDigi/interface/CSCStripDigiCollection.h"
00018 #include "DataFormats/FEDRawData/interface/FEDRawData.h"
00019 #include "DataFormats/FEDRawData/interface/FEDNumbering.h"
00020 #include "DataFormats/FEDRawData/interface/FEDRawDataCollection.h"
00021 #include "IORawData/CSCCommissioning/src/FileReaderDDU.h"
00022 #include "EventFilter/CSCRawToDigi/interface/CSCDDUEventData.h"
00023 #include "EventFilter/CSCRawToDigi/interface/CSCDCCEventData.h"
00024 #include "EventFilter/CSCRawToDigi/interface/CSCEventData.h"
00025 #include "EventFilter/CSCRawToDigi/interface/CSCDMBHeader.h"
00026 #include "OnlineDB/CSCCondDB/interface/CSCOldCrossTalkAnalyzer.h"
00027 #include "OnlineDB/CSCCondDB/interface/CSCxTalk.h"
00028 #include "CondFormats/CSCObjects/interface/CSCcrosstalk.h"
00029 #include "CondFormats/CSCObjects/interface/CSCPedestals.h"
00030 #include "OnlineDB/CSCCondDB/interface/CSCMap.h"
00031
00032 CSCOldCrossTalkAnalyzer::CSCOldCrossTalkAnalyzer(edm::ParameterSet const& conf) {
00033
00034 debug = conf.getUntrackedParameter<bool>("debug",false);
00035 eventNumber=0, Nddu=0,chamber=0;
00036 strip=0,misMatch=0,max1 =-9999999.,max2=-9999999.,adcMAX=-99999, min1=9999999.;
00037 layer=0,reportedChambers=0;
00038 length=1,myevt=0,flagRMS=-9,flagNoise=-9;
00039 aPeak=0.0,sumFive=0.0,maxPed=-99999.0,maxRMS=-99999.0;
00040 maxPeakTime=-9999.0, maxPeakADC=-999.0;
00041 minPeakTime=9999.0, minPeakADC=999.0;
00042 pedMean=0.0,evt=0,NChambers=0,myIndex=0;
00043
00044
00045 xtime = TH1F("time", "time", 50, 0, 500 );
00046 pulseshape_ch0=TH2F("Pulse Ch=0","All strips", 96,-100,500,100,-100,1100);
00047 pulseshape_ch1=TH2F("Pulse Ch=1","All strips", 96,-100,500,100,-100,1100);
00048 pulseshape_ch2=TH2F("Pulse Ch=2","All strips", 96,-100,500,100,-100,1100);
00049 pulseshape_ch3=TH2F("Pulse Ch=3","All strips", 96,-100,500,100,-100,1100);
00050 pulseshape_ch4=TH2F("Pulse Ch=4","All strips", 96,-100,500,100,-100,1100);
00051 pulseshape_ch5=TH2F("Pulse Ch=5","All strips", 96,-100,500,100,-100,1100);
00052 pulseshape_ch6=TH2F("Pulse Ch=6","All strips", 96,-100,500,100,-100,1100);
00053 pulseshape_ch7=TH2F("Pulse Ch=7","All strips", 96,-100,500,100,-100,1100);
00054 pulseshape_ch8=TH2F("Pulse Ch=8","All strips", 96,-100,500,100,-100,1100);
00055 pulseshape_ch9=TH2F("Pulse Ch=9","All strips", 96,-100,500,100,-100,1100);
00056 pulseshape_ch10=TH2F("Pulse Ch=10","All strips", 96,-100,500,100,-100,1100);
00057 pulseshape_ch11=TH2F("Pulse Ch=11","All strips", 96,-100,500,100,-100,1100);
00058 pulseshape_ch12=TH2F("Pulse Ch=12","All strips", 96,-100,500,100,-100,1100);
00059 ped_mean_all = TH1F("pedMean","Mean baseline noise", 100,300,900);
00060 maxADC = TH1F("maxADC","Peak ADC", 100,800,1300);
00061
00062
00063 for (int i=0;i<TOTALSTRIPS_oldxt;i++){
00064 new_xtalk_intercept_right[i] = -999.;
00065 new_xtalk_intercept_left[i] = -999.;
00066 new_xtalk_slope_right[i] = -999.;
00067 new_xtalk_slope_left[i] = -999.;
00068 new_rchi2[i] = -999.;
00069 new_lchi2[i] = -999.;
00070 newPeakTime[i] = -999.;
00071 newMeanPeakTime[i] = -999.;
00072 newPed[i] = 0 ;
00073 newRMS[i] = 0.0;
00074 newPeakRMS[i] = 0.0;
00075 newPeak[i] = 0.0;
00076 newPeakMin[i] = 0.0;
00077 newSumFive[i] = 0.0;
00078 }
00079
00080 for (int l=0; l<TIMEBINS_oldxt; l++){
00081 myTime[l] = 0.0;
00082 myADC[l] = 0.0;
00083 myTbin[l] = 0;
00084 }
00085
00086 for (int i=0;i<CHAMBERS_oldxt;i++){
00087 size[i] = 0;
00088 }
00089
00090 for (int iii=0;iii<DDU_oldxt;iii++){
00091 for (int i=0; i<CHAMBERS_oldxt; i++){
00092 for (int j=0; j<LAYERS_oldxt; j++){
00093 for (int k=0; k<STRIPS_oldxt; k++){
00094 for (int l=0; l<TIMEBINS_oldxt*PULSES_oldxt; l++){
00095 thetime[iii][i][j][k][l] = 0.0;
00096 thebins[iii][i][j][k][l] = 0 ;
00097 theadccountsc[iii][i][j][k][l] = 0 ;
00098 theadccountsl[iii][i][j][k][l] = 0 ;
00099 theadccountsr[iii][i][j][k][l] = 0 ;
00100 arrayOfPed[iii][i][j][k] = 0.;
00101 arrayOfPedSquare[iii][i][j][k] = 0.;
00102 arrayPed[iii][i][j][k] = 0.;
00103 arrayPeak[iii][i][j][k] = 0.;
00104 arrayPeakMin[iii][i][j][k] = 0.;
00105 arrayOfPeak[iii][i][j][k] = 0.;
00106 arrayOfPeakSquare[iii][i][j][k]= 0.;
00107 arraySumFive[iii][i][j][k] = 0.;
00108
00109 }
00110 }
00111 }
00112 }
00113 }
00114
00115
00116 for (int iii=0;iii<DDU_oldxt;iii++){
00117 for (int i=0; i<CHAMBERS_oldxt; i++){
00118 for (int j=0; j<LAYERS_oldxt; j++){
00119 for (int k=0; k<STRIPS_oldxt; k++){
00120 xtalk_intercept_left[iii][i][j][k] = -999.;
00121 xtalk_intercept_right[iii][i][j][k] = -999.;
00122 xtalk_slope_left[iii][i][j][k] = -999.;
00123 xtalk_slope_right[iii][i][j][k] = -999.;
00124 xtalk_chi2_left[iii][i][j][k] = -999.;
00125 xtalk_chi2_right[iii][i][j][k] = -999.;
00126 myPeakTime[iii][i][j][k] = 0.0 ;
00127 myMeanPeakTime[iii][i][j][k] = 0.0 ;
00128 array_meanPeakTime[iii][i][j][k] = -999.;
00129 }
00130 }
00131 }
00132 }
00133 }
00134
00135 void CSCOldCrossTalkAnalyzer::analyze(edm::Event const& e, edm::EventSetup const& iSetup) {
00136
00137 edm::Handle<CSCStripDigiCollection> strips;
00138 e.getByLabel("cscunpacker","MuonCSCStripDigi",strips);
00139 edm::Handle<FEDRawDataCollection> rawdata;
00140 e.getByType(rawdata);
00141 myevt=e.id().event();
00142
00143 for (int id=FEDNumbering::getCSCFEDIds().first;
00144 id<=FEDNumbering::getCSCFEDIds().second; ++id){
00145
00146
00148 const FEDRawData& fedData = rawdata->FEDData(id);
00149 if (fedData.size()){
00150
00152 CSCDCCEventData dccData((short unsigned int *) fedData.data());
00153
00154 const std::vector<CSCDDUEventData> & dduData = dccData.dduData();
00155
00156 evt++;
00157
00158 for (unsigned int iDDU=0; iDDU<dduData.size(); ++iDDU) {
00159
00161 const std::vector<CSCEventData> & cscData = dduData[iDDU].cscData();
00162 Nddu = dduData.size();
00163 reportedChambers += dduData[iDDU].header().ncsc();
00164 NChambers = cscData.size();
00165 int repChambers = dduData[iDDU].header().ncsc();
00166 std::cout << " Reported Chambers = " << repChambers <<" "<<NChambers<< std::endl;
00167 if (NChambers!=repChambers) { std::cout<< "misMatched size!!!" << std::endl; misMatch++; continue;}
00168
00169 for (int chamber = 0; chamber < NChambers; chamber++){
00170
00171 for (int layer = 1; layer <= 6; layer++){
00172
00173 std::vector<CSCStripDigi> digis = cscData[chamber].stripDigis(layer) ;
00174 const CSCDMBHeader * thisDMBheader = cscData[chamber].dmbHeader();
00175
00176 if (thisDMBheader->cfebAvailable()){
00177 dmbID[chamber] = cscData[chamber].dmbHeader()->dmbID();
00178 crateID[chamber] = cscData[chamber].dmbHeader()->crateID();
00179 if(crateID[chamber] == 255) continue;
00180
00181 for (unsigned int i=0; i<digis.size(); i++){
00182 size[chamber] = digis.size();
00183 int strip = digis[i].getStrip();
00184 std::vector<int> adc = digis[i].getADCCounts();
00185 if (adc[0]>0 && adc[1]>0) pedMean1 =(adc[0]+adc[1])/2;
00186 int offset = (evt-1)/PULSES_oldxt;
00187 int smain[5],splus[5],sminus[5];
00188 for(int s=0;s<5;s++) smain[s] = s*16+offset;
00189 for(int s=0;s<5;s++) splus[s] = s*16+offset+1;
00190 for(int s=0;s<5;s++) sminus[s] = s*16+offset-1;
00191 int iuse=-99;
00192 for(int s=0; s<5; s++) {if(strip-1==smain[s]) iuse=smain[s];}
00193 for(int s=0; s<5; s++) {if(strip-1==splus[s]) iuse=smain[s];}
00194 for(int s=0; s<5; s++) {if(strip-1==sminus[s]) iuse=smain[s];}
00195
00196 if(iuse!=-99){
00197
00198 for(unsigned int k=0;k<adc.size();k++){
00199
00200 time = (50. * k)-(((evt-1)%PULSES_oldxt)* 6.25)+116.5;
00201 pedMean =(adc[0]+adc[1])/2;
00202 ped_mean_all.Fill(pedMean);
00203 xtime.Fill(time);
00204
00205 if(chamber==0) pulseshape_ch0.Fill(time,adc[k]-pedMean);
00206 if(chamber==1) pulseshape_ch1.Fill(time,adc[k]-pedMean);
00207 if(chamber==2) pulseshape_ch2.Fill(time,adc[k]-pedMean);
00208 if(chamber==3) pulseshape_ch3.Fill(time,adc[k]-pedMean);
00209 if(chamber==4) pulseshape_ch4.Fill(time,adc[k]-pedMean);
00210 if(chamber==5) pulseshape_ch5.Fill(time,adc[k]-pedMean);
00211 if(chamber==6) pulseshape_ch6.Fill(time,adc[k]-pedMean);
00212 if(chamber==7) pulseshape_ch7.Fill(time,adc[k]-pedMean);
00213 if(chamber==8) pulseshape_ch8.Fill(time,adc[k]-pedMean);
00214 if(chamber==9) pulseshape_ch9.Fill(time,adc[k]-pedMean);
00215 if(chamber==10) pulseshape_ch10.Fill(time,adc[k]-pedMean);
00216 if(chamber==11) pulseshape_ch11.Fill(time,adc[k]-pedMean);
00217 if(chamber==12) pulseshape_ch12.Fill(time,adc[k]-pedMean);
00218
00219 myTime[k]=time;
00220 myADC[k]=adc[k];
00221 myTbin[k]=k;
00222
00223 aPeak = adc[3];
00224 if (max1 < aPeak) {
00225 max1 = aPeak;
00226 }
00227 sumFive = adc[2]+adc[3]+adc[4];
00228
00229 if (max2<sumFive){
00230 max2=sumFive;
00231 }
00232
00233 if (min1 > aPeak){
00234 min1=aPeak;
00235 }
00236
00237 if (pedMean1>0) maxADC.Fill(max1-pedMean1);
00238
00239 int kk=8*k-(evt-1)%PULSES_oldxt+19;
00240
00241 thebins[iDDU][chamber][layer-1][strip-1][kk] = 8*k-(evt-1)%PULSES_oldxt+19;
00242 thetime[iDDU][chamber][layer-1][strip-1][kk] = time;
00243
00244 if(iuse==strip-1) theadccountsc[iDDU][chamber][layer-1][iuse][kk] = adc[k];
00245 if(iuse==strip) theadccountsr[iDDU][chamber][layer-1][iuse][kk] = adc[k];
00246 if(iuse==strip-2) theadccountsl[iDDU][chamber][layer-1][iuse][kk] = adc[k];
00247 }
00248 }
00249
00250 arrayPed[iDDU][chamber][layer-1][strip-1] = pedMean1;
00251 arrayOfPed[iDDU][chamber][layer-1][strip-1] += pedMean1/evt;
00252 arrayOfPedSquare[iDDU][chamber][layer-1][strip-1] += pedMean1*pedMean1/(evt*evt) ;
00253 arrayPeak[iDDU][chamber][layer-1][strip-1] = max1-pedMean1;
00254 arrayPeakMin[iDDU][chamber][layer-1][strip-1] = min1;
00255 arrayOfPeak[iDDU][chamber][layer-1][strip-1] += max1-pedMean1;
00256 arrayOfPeakSquare[iDDU][chamber][layer-1][strip-1] += (max1-pedMean1)*(max1-pedMean1)/(evt*evt);
00257 arraySumFive[iDDU][chamber][layer-1][strip-1] = (max2-pedMean1)/(max1-pedMean1);
00258
00259 }
00260 }
00261 }
00262 }
00263
00264 eventNumber++;
00265 edm::LogInfo ("CSCOldCrossTalkAnalyzer") << "end of event number " << eventNumber;
00266
00267 }
00268 }
00269 }
00270 }
00271
00272 CSCOldCrossTalkAnalyzer::~CSCOldCrossTalkAnalyzer(){
00273 std::cout << "entering destructor " << std::endl;
00274
00275 Conv binsConv;
00276
00277 filein.open("../test/CSColdxtalk.cfg");
00278 filein.ignore(1000,'\n');
00279
00280 while(filein != NULL){
00281 lines++;
00282 getline(filein,PSet);
00283
00284 if (lines==2){
00285 name=PSet;
00286 }
00287 }
00288
00289 std::cout << "getting file name " << std::endl;
00290
00291
00292 std::string::size_type runNameStart = name.find("\"",0);
00293 std::string::size_type runNameEnd = name.find("raw",0);
00294 std::string::size_type rootStart = name.find("CrossTalk",0);
00295
00296 int nameSize = runNameEnd+2-runNameStart;
00297 int myRootSize = rootStart-runNameStart+8;
00298 std::string myname= name.substr(runNameStart+1,nameSize);
00299 std::string myRootName= name.substr(runNameStart+1,myRootSize);
00300 std::string myRootEnd = ".root";
00301 std::string myASCIIFileEnd = ".dat";
00302 std::string runFile= myRootName;
00303 std::string myRootFileName = runFile+myRootEnd;
00304 std::string myASCIIFileName= runFile+myASCIIFileEnd;
00305 const char *myNewName=myRootFileName.c_str();
00306 const char *myFileName=myASCIIFileName.c_str();
00307
00308 struct tm* clock;
00309 struct stat attrib;
00310 stat(myname.c_str(), &attrib);
00311 clock = localtime(&(attrib.st_mtime));
00312 std::string myTime=asctime(clock);
00313 std::ofstream myXtalkfile(myFileName,std::ios::app);
00314
00315
00316 cscmap *map = new cscmap();
00317
00318 TCalibOldCrossTalkEvt calib_evt;
00319 TFile calibfile(myNewName, "RECREATE");
00320 TTree calibtree("Calibration","Crosstalk");
00321 calibtree.Branch("EVENT", &calib_evt, "xtalk_slope_left/F:xtalk_slope_right/F:xtalk_int_left/F:xtalk_int_right/F:xtalk_chi2_left/F:xtalk_chi2_right/F:peakTime/F:strip/I:layer/I:cham/I:ddu/I:pedMean/F:pedRMS/F:peakRMS/F:maxADC/F:sum/F:id/I:flagRMS/I:flagNoise/I:MaxPed[13]/F:MaxRMS[13]/F:MaxPeakTime[13]/F:MinPeakTime[13]/F:MaxPeakADC[13]/F:MinPeakADC[13]/F");
00322 xtime.Write();
00323 ped_mean_all.Write();
00324 maxADC.Write();
00325 pulseshape_ch0.Write();
00326 pulseshape_ch1.Write();
00327 pulseshape_ch2.Write();
00328 pulseshape_ch3.Write();
00329 pulseshape_ch4.Write();
00330 pulseshape_ch5.Write();
00331 pulseshape_ch6.Write();
00332 pulseshape_ch7.Write();
00333 pulseshape_ch8.Write();
00334 pulseshape_ch9.Write();
00335 pulseshape_ch10.Write();
00336 pulseshape_ch11.Write();
00337 pulseshape_ch12.Write();
00338
00340
00341 float adc_ped_sub_left = -999.;
00342 float adc_ped_sub = -999.;
00343 float adc_ped_sub_right = -999.;
00344 int thebin;
00345 float sum=0.0;
00346 float mean=0;
00347
00348 for (int iii=0; iii<Nddu; iii++){
00349
00350 for (int i=0; i<NChambers; i++){
00351
00352
00353 int new_crateID = crateID[i];
00354 int new_dmbID = dmbID[i];
00355 int counter=0;
00356 std::cout<<" Crate: "<<new_crateID<<" and DMB: "<<new_dmbID<<std::endl;
00357 map->crate_chamber(new_crateID,new_dmbID,&chamber_id,&chamber_num,§or,&first_strip_index,&strips_per_layer,&chamber_index);
00358 std::cout<<"Data is for chamber:: "<< chamber_num<<" in sector: "<<sector<<std::endl;
00359
00360 calib_evt.id = chamber_num;
00361
00362 meanPedestal = 0.0;
00363 meanPeak = 0.0;
00364 meanPeakSquare=0.0;
00365 meanPedestalSquare = 0.;
00366 theRMS =0.0;
00367 thePedestal =0.0;
00368 theRSquare =0.0;
00369 thePeak =0.0;
00370 thePeakMin =0.0;
00371 thePeakRMS =0.0;
00372 theSumFive =0.0;
00373
00374 for (int j=0; j<LAYERS_oldxt; j++){
00375 mean=0.,sum=0.;
00376 for (int s=0; s<size[i]; s++) {
00377
00378 for (int m=0; m<3; m++){
00379 for (int n=0; n<120; n++){
00380 binsConv.convd[m][n] = 0.;
00381 binsConv.nconvd[m][n] = 0.;
00382 }
00383 }
00384
00385 for (int l=0; l < TIMEBINS_oldxt*PULSES_oldxt; l++){
00386 adc_ped_sub_left = theadccountsl[iii][i][j][s][l] - (theadccountsl[iii][i][j][s][0] + theadccountsl[iii][i][j][s][1])/2.;
00387 adc_ped_sub = theadccountsc[iii][i][j][s][l] - (theadccountsc[iii][i][j][s][0] + theadccountsc[iii][i][j][s][1])/2.;
00388 adc_ped_sub_right = theadccountsr[iii][i][j][s][l] - (theadccountsr[iii][i][j][s][0] + theadccountsr[iii][i][j][s][1])/2.;
00389 thebin=thebins[iii][i][j][s][l];
00390
00391 if (thebin >= 0 && thebin < 120){
00392 binsConv.convd[0][thebin] += adc_ped_sub_left;
00393 binsConv.nconvd[0][thebin] += 1.0;
00394
00395 binsConv.convd[1][thebin] += adc_ped_sub;
00396 binsConv.nconvd[1][thebin] += 1.0;
00397
00398 binsConv.convd[2][thebin] += adc_ped_sub_right;
00399 binsConv.nconvd[2][thebin] += 1.0;
00400
00401 }
00402 }
00403
00404 for (int m=0; m<3; m++){
00405 for (int n=0; n<120; n++){
00406 if(binsConv.nconvd[m][n]>1.0 && binsConv.nconvd[m][n] !=0.){
00407 binsConv.convd[m][n] = binsConv.convd[m][n]/binsConv.nconvd[m][n];
00408 }
00409 }
00410 }
00411
00412
00413 float xl_temp_a = 0.0;
00414 float xl_temp_b = 0.0;
00415 float minl_temp = 0.0;
00416 float xr_temp_a = 0.0;
00417 float xr_temp_b = 0.0;
00418 float minr_temp = 0.0;
00419 float pTime = 0.0;
00420
00421 binsConv.mkbins(50.);
00422 binsConv.convolution(&xl_temp_a, &xl_temp_b, &minl_temp, &xr_temp_a, &xr_temp_b, &minr_temp, &pTime);
00423 myPeakTime[iii][i][j][s] = pTime;
00424 sum=sum+myPeakTime[iii][i][j][s];
00425 mean = sum/size[i];
00426 }
00427
00428 int layer_id=chamber_num+j+1;
00429 if(sector==-100)continue;
00430
00431 for (int k=0; k<size[i]; k++){
00432
00433 for (int m=0; m<3; m++){
00434 for (int n=0; n<120; n++){
00435 binsConv.convd[m][n] = 0.;
00436 binsConv.nconvd[m][n] = 0.;
00437 }
00438 }
00439
00440 for (int l=0; l < TIMEBINS_oldxt*PULSES_oldxt; l++){
00441 adc_ped_sub_left = theadccountsl[iii][i][j][k][l] - (theadccountsl[iii][i][j][k][0] + theadccountsl[iii][i][j][k][1])/2.;
00442 adc_ped_sub = theadccountsc[iii][i][j][k][l] - (theadccountsc[iii][i][j][k][0] + theadccountsc[iii][i][j][k][1])/2.;
00443 adc_ped_sub_right = theadccountsr[iii][i][j][k][l] - (theadccountsr[iii][i][j][k][0] + theadccountsr[iii][i][j][k][1])/2.;
00444
00445 thebin=thebins[iii][i][j][k][l];
00446
00447 if (thebin >= 0 && thebin < 120){
00448 binsConv.convd[0][thebin] += adc_ped_sub_left;
00449 binsConv.nconvd[0][thebin] += 1.0;
00450
00451 binsConv.convd[1][thebin] += adc_ped_sub;
00452 binsConv.nconvd[1][thebin] += 1.0;
00453
00454 binsConv.convd[2][thebin] += adc_ped_sub_right;
00455 binsConv.nconvd[2][thebin] += 1.0;
00456
00457 }
00458 }
00459
00460 for (int m=0; m<3; m++){
00461 for (int n=0; n<120; n++){
00462 if(binsConv.nconvd[m][n]>1.0 && binsConv.nconvd[m][n] !=0.){
00463 binsConv.convd[m][n] = binsConv.convd[m][n]/binsConv.nconvd[m][n];
00464 }
00465 }
00466 }
00468
00470 float xl_temp_a = 0.;
00471 float xl_temp_b = 0.;
00472 float minl_temp = 0.;
00473 float xr_temp_a = 0.;
00474 float xr_temp_b = 0.;
00475 float minr_temp = 0.;
00476 float pTime = 0.;
00477
00478 binsConv.mkbins(50.);
00479 binsConv.convolution(&xl_temp_a, &xl_temp_b, &minl_temp, &xr_temp_a, &xr_temp_b, &minr_temp, &pTime);
00480
00481 if (k==0){
00482 xtalk_intercept_left[iii][i][j][k] = 0.0;
00483 xtalk_slope_left[iii][i][j][k] = 0.0;
00484 xtalk_chi2_left[iii][i][j][k] = 0.0;
00485
00486 xtalk_slope_right[iii][i][j][k] = xl_temp_b;
00487 xtalk_intercept_right[iii][i][j][k] = xl_temp_a;
00488 xtalk_chi2_right[iii][i][j][k] = minl_temp;
00489 myPeakTime[iii][i][j][k] = pTime;
00490
00491 if (isnan(xtalk_intercept_left[iii][i][j][k])) xtalk_intercept_left[iii][i][j][k] = 0.0;
00492 if (isnan(xtalk_intercept_right[iii][i][j][k])) xtalk_intercept_right[iii][i][j][k] = 0.0;
00493 if (isnan(xtalk_slope_left[iii][i][j][k])) xtalk_slope_left[iii][i][j][k] = 0.0;
00494 if (isnan(xtalk_slope_right[iii][i][j][k])) xtalk_slope_right[iii][i][j][k] = 0.0;
00495 if (isinf(xtalk_intercept_left[iii][i][j][k])) xtalk_intercept_left[iii][i][j][k] = 0.0;
00496 if (isinf(xtalk_intercept_right[iii][i][j][k])) xtalk_intercept_right[iii][i][j][k] = 0.0;
00497 if (isinf(xtalk_slope_left[iii][i][j][k])) xtalk_slope_left[iii][i][j][k] = 0.0;
00498 if (isinf(xtalk_slope_right[iii][i][j][k])) xtalk_slope_right[iii][i][j][k] = 0.0;
00499 }else if(k==size[i]-1){
00500 xtalk_intercept_right[iii][i][j][k] = 0.0;
00501 xtalk_slope_right[iii][i][j][k] = 0.0;
00502 xtalk_chi2_right[iii][i][j][k] = 0.0;
00503
00504 xtalk_intercept_left[iii][i][j][k] = xr_temp_a;
00505 xtalk_slope_left[iii][i][j][k] = xr_temp_b;
00506 xtalk_chi2_left[iii][i][j][k] = minr_temp;
00507 myPeakTime[iii][i][j][k] = pTime;
00508
00509 if (isnan(xtalk_intercept_left[iii][i][j][k])) xtalk_intercept_left[iii][i][j][k] = 0.0;
00510 if (isnan(xtalk_intercept_right[iii][i][j][k])) xtalk_intercept_right[iii][i][j][k] = 0.0;
00511 if (isnan(xtalk_slope_left[iii][i][j][k])) xtalk_slope_left[iii][i][j][k] = 0.0;
00512 if (isnan(xtalk_slope_right[iii][i][j][k])) xtalk_slope_right[iii][i][j][k] = 0.0;
00513 if (isinf(xtalk_intercept_left[iii][i][j][k])) xtalk_intercept_left[iii][i][j][k] = 0.0;
00514 if (isinf(xtalk_intercept_right[iii][i][j][k])) xtalk_intercept_right[iii][i][j][k] = 0.0;
00515 if (isinf(xtalk_slope_left[iii][i][j][k])) xtalk_slope_left[iii][i][j][k] = 0.0;
00516 if (isinf(xtalk_slope_right[iii][i][j][k])) xtalk_slope_right[iii][i][j][k] = 0.0;
00517 }else{
00518 xtalk_intercept_left[iii][i][j][k] = xl_temp_a;
00519 xtalk_intercept_right[iii][i][j][k] = xr_temp_a;
00520 xtalk_slope_left[iii][i][j][k] = xl_temp_b;
00521 xtalk_slope_right[iii][i][j][k] = xr_temp_b;
00522 xtalk_chi2_left[iii][i][j][k] = minl_temp;
00523 xtalk_chi2_right[iii][i][j][k] = minr_temp;
00524
00525 if (isnan(xtalk_intercept_left[iii][i][j][k])) xtalk_intercept_left[iii][i][j][k] = 0.0;
00526 if (isnan(xtalk_intercept_right[iii][i][j][k])) xtalk_intercept_right[iii][i][j][k] = 0.0;
00527 if (isnan(xtalk_slope_left[iii][i][j][k])) xtalk_slope_left[iii][i][j][k] = 0.0;
00528 if (isnan(xtalk_slope_right[iii][i][j][k])) xtalk_slope_right[iii][i][j][k] = 0.0;
00529 if (isinf(xtalk_intercept_left[iii][i][j][k])) xtalk_intercept_left[iii][i][j][k] = 0.0;
00530 if (isinf(xtalk_intercept_right[iii][i][j][k])) xtalk_intercept_right[iii][i][j][k] = 0.0;
00531 if (isinf(xtalk_slope_left[iii][i][j][k])) xtalk_slope_left[iii][i][j][k] = 0.0;
00532 if (isinf(xtalk_slope_right[iii][i][j][k])) xtalk_slope_right[iii][i][j][k] = 0.0;
00533
00534 myPeakTime[iii][i][j][k] = pTime;
00535 }
00536
00537 fff = (j*size[i])+k;
00538 float the_xtalk_left_a = xtalk_intercept_left[iii][i][j][k];
00539 float the_xtalk_right_a = xtalk_intercept_right[iii][i][j][k];
00540 float the_xtalk_left_b = xtalk_slope_left[iii][i][j][k];
00541 float the_xtalk_right_b = xtalk_slope_right[iii][i][j][k];
00542 float the_chi2_right = xtalk_chi2_right[iii][i][j][k];
00543 float the_chi2_left = xtalk_chi2_left[iii][i][j][k];
00544 float the_peakTime = myPeakTime[iii][i][j][k];
00545
00546 new_xtalk_intercept_right[fff] = the_xtalk_right_a ;
00547 new_xtalk_intercept_left[fff] = the_xtalk_left_a ;
00548 new_xtalk_slope_right[fff] = the_xtalk_right_b ;
00549 new_xtalk_slope_left[fff] = the_xtalk_left_b ;
00550 new_rchi2[fff] = the_chi2_right;
00551 new_lchi2[fff] = the_chi2_left;
00552 newPeakTime[fff] = the_peakTime;
00553 newMeanPeakTime[fff] = the_peakTime-mean;
00554
00555
00556 thePedestal = arrayPed[iii][i][j][k];
00557 meanPedestal = arrayOfPed[iii][i][j][k]/evt;
00558 newPed[fff] = meanPedestal;
00559 meanPedestalSquare = arrayOfPedSquare[iii][i][j][k]/evt;
00560 theRMS = sqrt(fabs(meanPedestalSquare - meanPedestal*meanPedestal));
00561
00562 newRMS[fff] = theRMS;
00563 theRSquare = (thePedestal-meanPedestal)*(thePedestal-meanPedestal)/(theRMS*theRMS*theRMS*theRMS);
00564 thePeak = arrayPeak[iii][i][j][k];
00565 thePeakMin = arrayPeakMin[iii][i][j][k];
00566 meanPeak = arrayOfPeak[iii][i][j][k]/evt;
00567 meanPeakSquare = arrayOfPeakSquare[iii][i][j][k]/evt;
00568 thePeakRMS = sqrt(fabs(meanPeakSquare - meanPeak*meanPeak));
00569 newPeakRMS[fff] = thePeakRMS;
00570 newPeak[fff] = thePeak;
00571 newPeakMin[fff] = thePeakMin;
00572
00573 theSumFive = arraySumFive[iii][i][j][k];
00574 newSumFive[fff]=theSumFive;
00575
00576 if(maxPed<thePedestal){
00577 maxPed=thePedestal;
00578 }
00579
00580 if(maxRMS<theRMS){
00581 maxRMS=theRMS;
00582 }
00583
00584 if(maxPeakTime<the_peakTime){
00585 maxPeakTime=the_peakTime;
00586 }
00587 if(minPeakTime>the_peakTime){
00588 minPeakTime=the_peakTime;
00589 }
00590
00591
00592 if (theRMS>1.5 && theRMS <6.0) flagRMS = 1;
00593 if (theRMS>6.0 && theRMS<20.0) flagRMS = 2;
00594 if (theRMS<1.5) flagRMS = 3;
00595 if (theRMS>=20.0) flagRMS = 4;
00596
00597 if (meanPedestal <50.) flagNoise = 2;
00598 if (meanPedestal >50. && meanPedestal<200.) flagNoise = 3;
00599 if (meanPedestal >1000. && meanPedestal<3000.) flagNoise = 4;
00600 if (meanPedestal >3000.) flagNoise = 5;
00601 if (meanPedestal >200. && meanPedestal<1000.) flagNoise = 1;
00602
00603
00604 counter++;
00605 myIndex = first_strip_index+counter-1;
00606 if (counter>size[i]*LAYERS_oldxt) counter=0;
00607 myXtalkfile <<layer_id<<" "<<myIndex<<" "<<k<<" "<<the_xtalk_right_b<<" "<<the_xtalk_right_a<<" "<<the_chi2_right<<" "<<the_xtalk_left_b<<" "<<the_xtalk_left_a<<" "<<the_chi2_left<<std::endl;
00608
00609 calib_evt.xtalk_slope_left = xtalk_slope_left[iii][i][j][k];
00610 calib_evt.xtalk_slope_right = xtalk_slope_right[iii][i][j][k];
00611 calib_evt.xtalk_int_left = xtalk_intercept_left[iii][i][j][k];
00612 calib_evt.xtalk_int_right = xtalk_intercept_right[iii][i][j][k];
00613 calib_evt.xtalk_chi2_left = xtalk_chi2_left[iii][i][j][k];
00614 calib_evt.xtalk_chi2_right = xtalk_chi2_right[iii][i][j][k];
00615 calib_evt.peakTime = myPeakTime[iii][i][j][k];
00616 calib_evt.cham = i;
00617 calib_evt.ddu = iii;
00618 calib_evt.layer = j;
00619 calib_evt.strip = k;
00620 calib_evt.flagRMS = flagRMS;
00621 calib_evt.flagNoise = flagNoise;
00622 calib_evt.pedMean = newPed[fff];
00623 calib_evt.pedRMS = newRMS[fff];
00624 calib_evt.peakRMS = newPeakRMS[fff];
00625 calib_evt.maxADC = newPeak[fff];
00626 calib_evt.sum = newSumFive[fff];
00627 calib_evt.MaxPed[i] = maxPed;
00628 calib_evt.MaxRMS[i] = maxRMS;
00629 calib_evt.MaxPeakTime[i] = maxPeakTime;
00630 calib_evt.MinPeakTime[i] = minPeakTime;
00631 calib_evt.MaxPeakADC[i] = newPeak[fff];
00632 calib_evt.MinPeakADC[i] = newPeakMin[fff];
00633
00634 calibtree.Fill();
00635
00636 }
00637 }
00638 }
00639 }
00640
00641 calibfile.Write();
00642 calibfile.Close();
00643 }
1.5.4