Back to home page

sPhenix code displayed by LXR

 
 

    


File indexing completed on 2025-08-06 08:16:10

0001 #include <cstdlib>
0002 #include <iostream>
0003 #include <map>
0004 #include <string>
0005 #include <vector>
0006 
0007 #include "TChain.h"
0008 #include "TFile.h"
0009 #include "TTree.h"
0010 #include "TString.h"
0011 #include "TObjString.h"
0012 #include "TSystem.h"
0013 #include "TROOT.h"
0014 #include "TTreeReader.h"
0015 #include "TTreeReaderValue.h"
0016 #include "TTreeReaderArray.h"
0017 
0018 void TPC_Channel_QA_Backup() {
0019   gROOT->SetBatch(kTRUE);
0020   std::ofstream outdata;
0021   outdata.open("noisyChannels.txt");
0022   if( !outdata ) { // file couldn't be opened
0023       cerr << "Error: file could not be opened" << endl;
0024       exit(1);
0025   }
0026   for (int q = 0; q < 24; q++){
0027     string sectorNumber;
0028     if (q < 10) sectorNumber = "0"+std::to_string(q);
0029     else sectorNumber = std::to_string(q);   
0030 
0031     string runNumber = "beam-00026540";
0032     
0033       string fileName = "/sphenix/user/jamesj3j3/tpc/sPHENIXProjects/pedestal-run-26540/TPC_ebdc"+sectorNumber+"_"+runNumber+"-0000.prdf_TPCRawDataTree_skip100.root";
0034     //    string fileName = "/sphenix/user/jamesj3j3/tpc/sPHENIXProjects/pedestal-run-26540/TPC_ebdc02_"+runNumber+"-0000.prdf_TPCRawDataTree_skip100.root"; 
0035 
0036 //    string mod[26] = {"R2","R2","R1","R1","R1","R3","R3","R3","R3","R3","R3","R2","R2","R1","R2","R2","R1","R1","R2","R2","R3","R3","R3","R3","R3","R3"}; 
0037 //    string slot[26] = {"5","6","1","3","2","12","10","11","9","8","7","1","2","4","8","7","6","5","4","3","1","3","2","4","6","5"}; 
0038     int mod_arr[26]={2,2,1,1,1,3,3,3,3,3,3,2,2,1,2,2,1,1,2,2,3,3,3,3,3,3};
0039     int slot_arr[26] = {5,6,1,3,2,12,10,11,9,8,7,1,2,4,8,7,6,5,4,3,1,3,2,4,6,5}; 
0040 
0041     TFile *truthInput = TFile::Open(fileName.c_str());
0042     TTreeReader myReader("SampleTree", truthInput);
0043 
0044     TTreeReaderValue<Int_t> nSamples_in(myReader, "nSamples");
0045     TTreeReaderValue<Int_t> fee_in(myReader, "fee");
0046     TTreeReaderArray<UShort_t> adcSamples_in(myReader, "adcSamples");
0047     TTreeReaderValue<Int_t> Channel_in(myReader, "Channel");
0048 
0049 //Adding the output tree
0050     TString* outputfilename=new TString("/sphenix/user/jamesj3j3/tpc/sPHENIXProjects/pedestal-run-26540/outputfile_TPC_ebdc"+sectorNumber+"_"+runNumber+"_skip100.root");
0051     TFile* outputfile=new TFile(outputfilename->Data(),"recreate");
0052     TTree* outputTree=new TTree("outputTree","outputTree");
0053     Int_t isAlive=1.0; // 1 if the channel is working properly, 0 if no signal(number of samples is 0, adc value is 0 or nan), pedestal above 200 or below 10
0054     Float_t pedMean; // average pedestal value over all samples
0055     Float_t pedStdi; // pedestal standard deviation over all samples
0056     Int_t chan; // channel number
0057     Int_t fee; // fee number
0058     Int_t module; // module number (e.g. R1, R2, R3)
0059     Int_t slot; // slot number
0060     outputTree->Branch("isAlive",&isAlive,"isAlive/I"); 
0061     outputTree->Branch("pedMean",&pedMean,"pedMean/F");
0062         outputTree->Branch("pedStdi",&pedStdi,"pedStdi/F");
0063         outputTree->Branch("chan",&chan,"chan/I");
0064         outputTree->Branch("fee",&fee,"fee/I");
0065         outputTree->Branch("module",&module,"module/I");
0066     outputTree->Branch("slot",&slot,"slot/I");
0067 
0068     
0069 
0070     string dChan = "Dead Channels (Sector "+sectorNumber+", Run "+runNumber+");Channel Number + 256*FEE Number;Instances"; 
0071 
0072     //TH1F *ave_adc = new TH1F("ave_adc","Channel Pedestals;Channel Number + 256*FEE Number;Pedestal ADC",256*26,-0.5,255.5+256*25);
0073     TH1F *dead_channel = new TH1F("dead_channel",dChan.c_str(),256*26,-0.5,255.5+256*25);
0074     //TH1F *rms_adc = new TH1F("rms_adc","Pedestal RMS;Channel Number + 256*FEE Number;Pedestal RMS",256*26,-0.5,255.5+256*25);
0075   
0076     std::vector<float> channels;
0077     std::vector<float> pedestal;
0078     std::vector<float> std_data;
0079 
0080     Float_t ave_adc_fee_channel[26][256];
0081     Float_t std_adc_fee_channel[26][256];
0082     Float_t counts_adc_fee_channel[26][256];
0083     Int_t alive_array_fee_channel[26][256];
0084     for(int fee_no=0;fee_no<26;fee_no++){
0085         for(int channel_no=0;channel_no<256;channel_no++)
0086     {
0087             ave_adc_fee_channel[fee_no][channel_no]=0.0;
0088             std_adc_fee_channel[fee_no][channel_no]=0.0;
0089             counts_adc_fee_channel[fee_no][channel_no]=0.0;
0090         alive_array_fee_channel[fee_no][channel_no]=1;
0091         }
0092     }
0093 
0094     //Int_t channels[256];
0095     while(myReader.Next())
0096         {
0097             //cout<<"fee="<<*fee<<"Channel="<<*Channel<<endl;
0098             if(*nSamples_in==0) 
0099         {
0100         dead_channel->Fill(*Channel_in+ 256*(*fee_in)); 
0101         alive_array_fee_channel[*fee_in][*Channel_in]=0; 
0102         continue;
0103         }
0104             
0105             bool dead = false;
0106             for(int adc_sam_no=0;adc_sam_no<*nSamples_in;adc_sam_no++){
0107                 //if (*fee_in == 25) channels[*Channel_in] += adcSamples_in[adc_sam_no];
0108                 if (adcSamples_in[adc_sam_no] == 0 || TMath::IsNaN(float(adcSamples_in[adc_sam_no]))) {
0109             dead_channel->Fill(*Channel_in+ 256*(*fee_in));dead=true;break;
0110             alive_array_fee_channel[*fee_in][*Channel_in]=0;
0111         }
0112             }
0113             if (dead) {continue;}
0114 
0115             for(int adc_sam_no=0;adc_sam_no<*nSamples_in;adc_sam_no++){
0116                 ave_adc_fee_channel[*fee_in][*Channel_in]+=adcSamples_in[adc_sam_no];
0117                 std_adc_fee_channel[*fee_in][*Channel_in]+=pow(adcSamples_in[adc_sam_no],2);
0118                 counts_adc_fee_channel[*fee_in][*Channel_in]+=1.0;
0119             }
0120 
0121         }
0122 
0123     for(int fee_no=0;fee_no<26;fee_no++){
0124           for(int channel_no=0;channel_no<256;channel_no++){
0125             if(counts_adc_fee_channel[fee_no][channel_no] != 0.0) 
0126             {
0127                 float temp1 = ave_adc_fee_channel[fee_no][channel_no]/counts_adc_fee_channel[fee_no][channel_no];
0128                 float temp2 = std_adc_fee_channel[fee_no][channel_no]/counts_adc_fee_channel[fee_no][channel_no];
0129                 ave_adc_fee_channel[fee_no][channel_no] = temp1;
0130                 std_adc_fee_channel[fee_no][channel_no] = temp2;
0131             }
0132             else
0133         {
0134         dead_channel->Fill(channel_no+ 256*(fee_no));
0135                 //std::cout << "FEE: " << fee_no << ", channel: " << channel_no << std::endl;
0136                 ave_adc_fee_channel[fee_no][channel_no] = 0.0;
0137                 std_adc_fee_channel[fee_no][channel_no] = 0.0;
0138         alive_array_fee_channel[fee_no][channel_no]=0;
0139             }
0140             
0141         if(sqrt(std_adc_fee_channel[fee_no][channel_no] - pow(ave_adc_fee_channel[fee_no][channel_no],2)) > 4)
0142         {
0143               outdata << "Sector " << q << ", FEE " << fee_no << ", Channel " << channel_no << std::endl;
0144         }
0145             if(ave_adc_fee_channel[fee_no][channel_no] > 200)
0146             {
0147         dead_channel->Fill(channel_no+ 256*(fee_no));
0148         std::cout << "(Large Pedestal pedestal>200) FEE: " << fee_no << ", channel: " << channel_no << ", " << ave_adc_fee_channel[fee_no][channel_no] << std::endl;
0149         alive_array_fee_channel[fee_no][channel_no]=0;
0150             }
0151             else if(ave_adc_fee_channel[fee_no][channel_no] < 10)
0152             {
0153         dead_channel->Fill(channel_no+ 256*(fee_no));
0154         std::cout << "(Small Pedestal pedestal<10) FEE: " << fee_no << ", channel: " << channel_no << ", " << ave_adc_fee_channel[fee_no][channel_no] << std::endl;
0155         alive_array_fee_channel[fee_no][channel_no]=0;
0156             }
0157             channels.push_back(channel_no+256.*(fee_no));
0158             pedestal.push_back(ave_adc_fee_channel[fee_no][channel_no]);
0159             std_data.push_back(sqrt(std_adc_fee_channel[fee_no][channel_no] - pow(ave_adc_fee_channel[fee_no][channel_no],2)));
0160             //ave_adc->SetBinContent(channel_no*(fee_no+1),ave_adc_fee_channel[fee_no][channel_no]);
0161             //std::cout << ave_adc_fee_channel[fee_no][channel_no] << std::endl;
0162 
0163         pedMean=ave_adc_fee_channel[fee_no][channel_no];
0164         pedStdi=sqrt(std_adc_fee_channel[fee_no][channel_no] - pow(ave_adc_fee_channel[fee_no][channel_no],2));
0165         isAlive=alive_array_fee_channel[fee_no][channel_no];
0166             chan = channel_no;
0167             //if(*nSamples_in==0 || pedMean_o>200 || pedMean_o<10){isAlive=false};
0168         fee=fee_no;
0169         module=mod_arr[fee_no];
0170         slot=slot_arr[fee_no];
0171         outputTree->Fill();
0172         
0173 
0174           }
0175       }
0176       //for (int l = 0; l < 256; l++)
0177       //{
0178       //    std::cout << "Channel " << l << ": " << ave_adc_fee_channel[25][l] << std::endl;
0179       //}
0180 
0181       TGraph *ave_adc = new TGraph(256*26,channels.data(),pedestal.data());
0182       TGraph *rms_adc = new TGraph(256*26,channels.data(),std_data.data());
0183       
0184       TCanvas *c1 = new TCanvas("c1","c1",2000,1200);
0185       c1->Divide(2,2);
0186       c1->cd(1);
0187       auto tl = new TLine(); tl->SetLineColor(kRed);
0188       dead_channel->Draw();
0189       for (int k = 0; k < 26; k++)
0190       {
0191         tl->DrawLine(256*(k+1),0,256*(k+1),dead_channel->GetMaximum());
0192         TLatex *pt = new TLatex(25+k*256,dead_channel->GetMaximum()*5/6,std::to_string(k).c_str());
0193         pt->SetTextSize(0.04);
0194     pt->Draw();
0195         //TLatex *pt2 = new TLatex(25+k*256,dead_channel->GetMaximum()*7.5/10,mod[k].c_str());
0196         //pt2->SetTextSize(0.03);
0197     //pt2->Draw();
0198         //TLatex *pt3 = new TLatex(25+k*256,dead_channel->GetMaximum()*4/6,slot[k].c_str());
0199         //pt3->SetTextSize(0.04);
0200     //pt3->Draw();
0201       }
0202       c1->cd(2);
0203       ave_adc->SetTitle("Mean Pedestal Value;Channel Number + 256*FEE Number;Mean Pedestal ADC");
0204       ave_adc->SetMarkerStyle(20);
0205       ave_adc->SetMarkerColor(kBlue);
0206       ave_adc->SetMarkerSize(0.25);
0207       ave_adc->GetXaxis()->SetRangeUser(0,256*26);
0208       ave_adc->GetYaxis()->SetRangeUser(0,150);
0209       ave_adc->Draw("AP");
0210       for (int k = 0; k < 26; k++)
0211       {
0212         tl->DrawLine(256*(k+1),0,256*(k+1),150);
0213         TLatex *pt = new TLatex(25+k*256,135,std::to_string(k).c_str());
0214         pt->SetTextSize(0.04);
0215     pt->Draw();
0216     /*
0217         //TLatex *pt2 = new TLatex(25+k*256,125,mod[k].c_str());
0218         //pt2->SetTextSize(0.03);
0219     //pt2->Draw();
0220         TLatex *pt3 = new TLatex(25+k*256,115,slot[k].c_str());
0221         pt3->SetTextSize(0.04);
0222     pt3->Draw();
0223       */
0224       }
0225       c1->cd(3);
0226       rms_adc->SetTitle("Pedestal Standard Deviation;Channel Number + 256*FEE Number;Pedestal Standard Deviation");
0227       rms_adc->SetMarkerStyle(20);
0228       rms_adc->SetMarkerColor(kBlue);
0229       rms_adc->SetMarkerSize(0.25);
0230       rms_adc->GetXaxis()->SetRangeUser(0,256*26);
0231       rms_adc->GetYaxis()->SetRangeUser(0,15);
0232       rms_adc->Draw("AP");
0233       for (int k = 0; k < 26; k++)
0234       {
0235         tl->DrawLine(256*(k+1),0,256*(k+1),15);
0236         TLatex *pt = new TLatex(25+k*256,13.5,std::to_string(k).c_str());
0237         pt->SetTextSize(0.04);
0238     pt->Draw();
0239         //TLatex *pt2 = new TLatex(25+k*256,12.5,mod[k].c_str());
0240         //pt2->SetTextSize(0.03);
0241     //pt2->Draw();
0242         //TLatex *pt3 = new TLatex(25+k*256,11.5,slot[k].c_str());
0243         //pt3->SetTextSize(0.04);
0244     //pt3->Draw();
0245       }
0246  
0247       string saveName = "images/S"+sectorNumber+"-"+runNumber+".png";
0248       c1->SaveAs(saveName.c_str());
0249       c1->Close();
0250   outputfile->cd();
0251   outputTree->Write();
0252   outputfile->Close();
0253   }
0254   outdata.close();
0255 }