ARA ROOT v3.10 Software

calibration/ATRI/voltageCalib/getADCmVConversionTest.cxx

00001 
00002 
00003 
00004 
00005 
00006 //AraRoot Includes
00007 #include "UsefulAtriStationEvent.h"
00008 #include "RawAtriStationEvent.h"
00009 #include "araSoft.h"
00010 
00011 
00012 //Root Includes
00013 #include "TTree.h"
00014 #include "TFile.h"
00015 #include "TH1.h"
00016 #include "TF1.h"
00017 #include "TTree.h"
00018 #include "TMath.h"
00019 #include "TCanvas.h"
00020 #include "TChain.h"
00021 #include "TGraph.h"
00022 #include "TMultiGraph.h"
00023 
00024 
00025 //Standard Includes
00026 #include <iostream>
00027 #include <fstream>
00028 
00029 
00030 #define MAX_ADC 4096
00031 
00032 //Prototype Functions
00033 Int_t getADCmVConversion(char *baseName, Int_t runLow, Int_t runHigh, Int_t thisBlock, Int_t thisDda);
00034 TGraph *getZeroSubGraph(TGraph *gr, Int_t factor);
00035 TGraph *addGraphs(TGraph *grOne, TGraph *grTwo);
00036 
00037 int main(int argc, char **argv)
00038 {
00039   Int_t runLow, runHigh, block, dda;
00040   char baseName[FILENAME_MAX];
00041 
00042   if(argc<6) {
00043     std::cerr << "Usage: " << argv[0] << " <baseDir> <runLow> <runHigh> <block> <dda> \n";
00044     return 1;
00045   }
00046   sprintf(baseName, argv[1]);
00047   runLow=atoi(argv[2]);
00048   runHigh=atoi(argv[3]);
00049   block=atoi(argv[4]);
00050   dda=atoi(argv[5]);
00051   getADCmVConversion(baseName, runLow, runHigh, block, dda);
00052 
00053 }
00054 
00055 
00056 
00057 
00058 
00059 
00060 Int_t getADCmVConversion(char *baseName, Int_t runLow, Int_t runHigh, Int_t thisBlock, Int_t thisDda){
00061 
00062   TChain *chain = new TChain("maxPosNegTree");
00063   char runName[FILENAME_MAX];
00064   for(int run=runLow; run<=runHigh;run++){
00065     sprintf(runName, "%s/root/run%i/adcSampleBlock.root", baseName, run);
00066     chain->Add(runName);
00067   }
00068   char outFileName[FILENAME_MAX];
00069   sprintf(outFileName, "%s/root/voltageCalibThursday/ADCmVConversionTest_block%idda%i.root", baseName, thisBlock, thisDda);
00070   
00071   TFile *outFile = new TFile(outFileName, "RECREATE");
00072 
00073   Int_t dda=0, chan=0, sample=0, block=0, posNeg=0, ndf=0;
00074   Double_t maxPosValue=0, maxNegValue=0, voltage=0;
00075   Double_t p0=0, p1=0, p2=0, p3=0, p4=0, chiSq=0;
00076   
00077 
00078   
00079   char exp[100];
00080   char cut[100];
00081   char grName[100];
00082   TGraph *grTemp;
00083   TCanvas *tempCan;
00084 
00085   TF1* funcPos = new TF1("myPosPol3", "pol3");
00086   TF1* funcNeg = new TF1("myNegPol1", "pol1");
00087   funcPos->SetParLimits(1, 0,2);
00088   funcPos->SetParLimits(2, -1,+1);
00089   funcPos->SetParLimits(3, -1,+1);
00090   funcNeg->SetParLimits(1, 0,+2);
00091   funcPos->FixParameter(0,0);
00092   funcNeg->FixParameter(0,0);
00093         
00094 
00095 
00096 
00097   TTree *fitTree = new TTree("fitTree", "Tree of fit parameters");
00098   fitTree->Branch("dda", &dda, "dda/I");
00099   fitTree->Branch("chan", &chan, "chan/I");
00100   fitTree->Branch("sample", &sample, "sample/I");
00101   fitTree->Branch("block", &block, "block/I");
00102   fitTree->Branch("ndf", &ndf, "ndf/I");
00103   fitTree->Branch("posNeg", &posNeg, "posNeg/I");
00104 
00105   fitTree->Branch("p0", &p0, "p0/D");
00106   fitTree->Branch("p1", &p1, "p1/D");
00107   fitTree->Branch("p2", &p2, "p2/D");
00108   fitTree->Branch("p3", &p3, "p3/D");
00109   fitTree->Branch("chiSq", &chiSq, "chiSq/D");
00110 
00111 
00112 
00113   block=thisBlock;
00114   dda=thisDda;
00115   //  for(dda=0;dda<1;dda++){
00116     for(chan=0;chan<1;chan++){
00117       //      fprintf(stderr, "dda %i chan %i ", dda, chan);
00118       for(sample=40;sample<45;sample++){
00119 
00120         if((dda==1 || dda==2)&&(chan>3)) continue;
00121         //Reset starting point of parameters
00122         funcPos->SetParameter(1, 1.0);
00123         funcPos->SetParameter(2, -0.001);
00124         funcPos->SetParameter(3, 0);
00125         funcNeg->SetParameter(1, 0.55);
00126 
00127         //      fprintf(stderr, "%i ", sample);
00128         sprintf(exp, "voltage/2:maxPosValue");
00129         sprintf(cut, "dda==%i&&chan==%i&&block==%i&&sample==%i", dda, chan, block, sample);
00130         sprintf(grName, "grPos_%i_%i_%i_%i", dda, chan, block, sample);
00131         tempCan = new TCanvas();
00132         chain->Draw(exp, cut);
00133         grTemp = (TGraph*)tempCan->GetPrimitive("Graph");
00134         //      TGraph *grPosTemp = (TGraph*) grTemp->Clone();
00135         TGraph *grPosTemp = getZeroSubGraph(grTemp, +1);
00136 
00137         grPosTemp->SetName(grName);
00138         grPosTemp->Fit("myPosPol3", "QBW");
00139         p0=funcPos->GetParameter(0);
00140         p1=funcPos->GetParameter(1);
00141         p2=funcPos->GetParameter(2);
00142         p3=funcPos->GetParameter(3);
00143         chiSq =funcPos->GetChisquare();
00144         ndf = funcPos->GetNDF();
00145         p4=0;
00146         posNeg=+1;
00147         printf("%i %i %i %i %f %f\n", dda, chan, block, sample, p1, p2);
00148         fitTree->Fill();
00149         delete tempCan;
00150         
00151         
00152         tempCan = new TCanvas();
00153         sprintf(exp, "-1*voltage/2:maxNegValue");
00154         sprintf(grName, "grNeg_%i_%i_%i_%i", dda, chan, block, sample);   
00155         chain->Draw(exp, cut);
00156         grTemp = (TGraph*)tempCan->GetPrimitive("Graph");
00157         //      TGraph *grNegTemp = (TGraph*) grTemp->Clone();
00158         TGraph *grNegTemp = getZeroSubGraph(grTemp, -1);
00159 
00160         grNegTemp->SetName(grName);
00161         grNegTemp->Fit("myNegPol1", "QBW");
00162         p0=funcNeg->GetParameter(0);
00163         p1=funcNeg->GetParameter(1);
00164         chiSq=funcNeg->GetChisquare();
00165         p2=0;
00166         p3=0;
00167         p4=0;
00168         // p0=grNegTemp->GetFunction("myNegPol1")->GetParameter(0);
00169         // p1=grNegTemp->GetFunction("myNegPol1")->GetParameter(1);
00170         // chiSq = grNegTemp->GetFunction("myNegPol1")->GetChisquare();
00171 
00172         posNeg=-1;
00173         printf("%i %i %i %i %f %f\n", dda, chan, block, sample, p1, p2);
00174         fitTree->Fill();
00175         delete tempCan;
00176         
00177 
00178         outFile->cd();
00179         grPosTemp->Write();
00180         grNegTemp->Write();
00181         //      grPosNeg->Write();
00182 
00183         //      delete grPosNeg;
00184         delete grNegTemp;       
00185         delete grPosTemp;
00186 
00187 
00188       }
00189       fprintf(stderr, "\n");
00190     }
00191     //  }
00192 
00193   outFile->cd();
00194   fitTree->Write();
00195 
00196   fprintf(stderr, "Done!\n");
00197 
00198   delete funcPos;
00199   delete funcNeg;
00200 
00201 
00202   return 0;
00203 
00204 }
00205 
00206 
00207 
00208 TGraph *getZeroSubGraph(TGraph *gr, Int_t factor){
00209   Int_t numEntries = gr->GetN();
00210   Double_t *yVals=gr->GetY();
00211   Double_t *xVals=gr->GetX();
00212   
00213   if(numEntries<=0) return NULL;
00214   Double_t *xValsNew = new Double_t[numEntries];
00215   Double_t xValZero=0;
00216 
00217   for(int entry=0;entry<numEntries;entry++){
00218     if(yVals[entry]==0){
00219       xValZero=xVals[entry];
00220       break;
00221     }
00222   }
00223   for(int entry=0;entry<numEntries;entry++){
00224     Double_t newX = TMath::Sqrt( TMath::Power(xVals[entry],2) - TMath::Power(xValZero,2) );
00225     xValsNew[entry] = factor*newX;
00226   }
00227 
00228   TGraph *grNew = new TGraph(numEntries, xValsNew, yVals);
00229 
00230 
00231   delete [] xValsNew;
00232 
00233   return grNew;
00234 
00235 }
00236 
00237 TGraph *addGraphs(TGraph *grOne, TGraph *grTwo){
00238   Int_t numEntriesOne = grOne->GetN();
00239   Int_t numEntriesTwo = grTwo->GetN();
00240   Double_t *yVals1 = grOne->GetY();
00241   Double_t *xVals1 = grOne->GetX();
00242   Double_t *yVals2 = grTwo->GetY();
00243   Double_t *xVals2 = grTwo->GetX();
00244 
00245   Double_t *yValsNew = new Double_t[numEntriesOne+numEntriesTwo];
00246   Double_t *xValsNew = new Double_t[numEntriesOne+numEntriesTwo];
00247 
00248   for(int entry = 0;entry<numEntriesOne;entry++){
00249     xValsNew[entry] = xVals1[entry];
00250     yValsNew[entry] = yVals1[entry];
00251   }
00252   for(int entry = 0;entry<numEntriesTwo;entry++){
00253     xValsNew[entry+numEntriesOne] = xVals2[entry];
00254     yValsNew[entry+numEntriesOne] = yVals2[entry];
00255   }
00256 
00257   TGraph *grNew = new TGraph(numEntriesOne+numEntriesTwo, xValsNew, yValsNew);
00258 
00259   delete [] yValsNew;
00260   delete [] xValsNew;
00261 
00262   return grNew;
00263 }

Generated on Tue Jul 16 16:58:02 2013 for ARA ROOT v3.10 Software by doxygen 1.4.7