analysis/deltaTPulses.cxx
00001 00002 00003 00004 00005 00006 00007 00008 00009 00010 //Includes 00011 #include <iostream> 00012 00013 //AraRoot Includes 00014 #include "RawAtriStationEvent.h" 00015 #include "UsefulAraStationEvent.h" 00016 #include "UsefulAtriStationEvent.h" 00017 #include "AraGeomTool.h" 00018 #include "AraStationInfo.h" 00019 #include "AraAntennaInfo.h" 00020 #include "AraCalAntennaInfo.h" 00021 00022 00023 #include "FFTtools.h" 00024 00025 //ROOT Includes 00026 #include "TTree.h" 00027 #include "TMath.h" 00028 #include "TFile.h" 00029 #include "TGraph.h" 00030 00031 00032 RawAtriStationEvent *rawAtriEvPtr; 00033 UsefulAtriStationEvent *realAtriEvPtr; 00034 Double_t interSample=0.1; 00035 Double_t getMaxX(TGraph *gr); 00036 00037 int main(int argc, char **argv) 00038 { 00039 00040 if(argc<3) { 00041 std::cout << "Usage\n" << argv[0] << " <input file> <output file>\n"; 00042 return 0; 00043 } 00044 00045 printf("------------------------------------------------------------------------------------------------------------------------------------------\n"); 00046 printf("\t%s\n", argv[0]); 00047 printf("------------------------------------------------------------------------------------------------------------------------------------------\n"); 00048 00049 printf("input file %s\n", argv[1]); 00050 printf("out file %s\n", argv[2]); 00051 00052 TFile *fpOut = new TFile(argv[2], "RECREATE"); 00053 00054 TTree *outTree = new TTree("tDeltaT", "tree of cross correlations"); 00055 Double_t deltaT=0; 00056 UInt_t thisTimeStamp=0, pps=0; 00057 Int_t chan1=0, chan2=0; 00058 AraAntennaInfo *antInfo1=0, *antInfo2=0; 00059 00060 outTree->Branch("deltaT", &deltaT, "deltaT/D"); 00061 outTree->Branch("timeStamp", &thisTimeStamp, "thisTimeStamp/I"); 00062 outTree->Branch("pps", &pps, "pps/I"); 00063 outTree->Branch("chan1", &chan1, "chan1/I"); 00064 outTree->Branch("chan2", &chan2, "chan2/I"); 00065 outTree->Branch("antInfo1","AraAntennaInfo",&antInfo1); 00066 outTree->Branch("antInfo2","AraAntennaInfo",&antInfo2); 00067 00068 00069 TFile *fp = TFile::Open(argv[1]); 00070 if(!fp) { 00071 std::cerr << "Can't open file\n"; 00072 return -1; 00073 } 00074 TTree *eventTree = (TTree*) fp->Get("eventTree"); 00075 if(!eventTree) { 00076 std::cerr << "Can't find eventTree\n"; 00077 return -1; 00078 } 00079 00080 eventTree->SetBranchAddress("event", &rawAtriEvPtr); 00081 Long64_t numEntries=eventTree->GetEntries(); 00082 Long64_t starEvery=numEntries/80; 00083 if(starEvery==0) starEvery++; 00084 00085 TGraph *grChan1; 00086 TGraph *grChan2; 00087 TGraph *grChan1Int; 00088 TGraph *grChan2Int; 00089 TGraph * grCorr; 00090 00091 00092 00093 eventTree->GetEntry(0); 00094 00095 Int_t stationId = rawAtriEvPtr->stationId; 00096 AraStationInfo *tempStationInfo = new AraStationInfo(stationId); 00097 00098 for(Long64_t event=0;event<numEntries;event++) { 00099 if(event%starEvery==0) { 00100 std::cerr << "*"; 00101 } 00102 00103 //This line gets the RawIcrr or RawAtri Event 00104 eventTree->GetEntry(event); 00105 00106 thisTimeStamp = rawAtriEvPtr->timeStamp; 00107 pps = rawAtriEvPtr->ppsNumber; 00108 00109 if(!(rawAtriEvPtr->isCalpulserEvent())) continue; 00110 00111 realAtriEvPtr = new UsefulAtriStationEvent(rawAtriEvPtr, AraCalType::kFirstCalibPlusCables); 00112 00113 for(chan1=0;chan1<ANTS_PER_ATRI;chan1++){ 00114 antInfo1 = tempStationInfo->getAntennaInfo(chan1); 00115 00116 for(chan2=chan1+1;chan2<ANTS_PER_ATRI;chan2++){ 00117 00118 antInfo2 = tempStationInfo->getAntennaInfo(chan2); 00119 00120 grChan1 = realAtriEvPtr->getGraphFromRFChan(chan1); 00121 grChan2 = realAtriEvPtr->getGraphFromRFChan(chan2); 00122 00123 grChan1Int=FFTtools::getInterpolatedGraph(grChan1,interSample); 00124 00125 grChan2Int=FFTtools::getInterpolatedGraph(grChan2,interSample); 00126 00127 grCorr = FFTtools::getInterpolatedCorrelationGraph(grChan1Int,grChan2Int,interSample); 00128 // grCorr=FFTtools::getCorrelationGraph(grChan1Int,grChan2Int); 00129 00130 00131 deltaT = getMaxX(grCorr); 00132 00133 outTree->Fill(); 00134 00135 delete grChan1Int; 00136 delete grChan2Int; 00137 delete grCorr; 00138 delete grChan1; 00139 delete grChan2; 00140 } 00141 } 00142 delete realAtriEvPtr; 00143 00144 } 00145 std::cerr << "\n"; 00146 00147 fpOut->Write(); 00148 00149 00150 } 00151 00152 00153 Double_t getMaxX(TGraph *gr){ 00154 Int_t size = gr->GetN(); 00155 Double_t *xVals = gr->GetX(); 00156 Double_t *yVals = gr->GetY(); 00157 00158 Double_t maxX=0; 00159 Double_t maxY=0; 00160 for(int n=0;n<size;n++){ 00161 if(TMath::Abs(yVals[n])>TMath::Abs(maxY)){ 00162 maxX=xVals[n]; 00163 maxY=yVals[n]; 00164 00165 } 00166 00167 00168 } 00169 00170 return maxX; 00171 }
Generated on Tue Jul 16 16:58:01 2013 for ARA ROOT v3.10 Software by
1.4.7