import ROOT, glob, math, sys, os
from ROOT import std

maxDist = 3

#argument parsing:
for iarg in range(len(sys.argv)) :
    arg = sys.argv[iarg]
    if arg == "--prelim" :
        preliminary = 1
    if arg == "-f" :
        tfile = sys.argv[iarg+1]
    if arg == "-d" :
        maxDist = float (sys.argv[iarg+1])
#finished argument parsing


ROOT.gROOT.ProcessLine(".x rootlogon.C") 
ROOT.gROOT.ProcessLine('gStyle->SetPalette(1)')
ROOT.gROOT.ProcessLine('gROOT->SetStyle("Plain")')
ROOT.gROOT.ProcessLine('gStyle->SetOptStat(0)')

#histograms:
resX = ROOT.TH1D("resX","Laser matching residual X", 200,-maxDist, maxDist)
resX.GetXaxis().SetTitle("cm")
resY = ROOT.TH1D("resY","Laser matching residual Y", 200,-maxDist, maxDist)
resY.GetXaxis().SetTitle("cm")
resZ = ROOT.TH1D("resZ","Laser matching residual Z", 200,-maxDist, maxDist)
resZ.GetXaxis().SetTitle("cm")

reslX = ROOT.TH1D("reslX","Laser track residual X", 200,-maxDist, maxDist)
resX.GetXaxis().SetTitle("cm")
reslY = ROOT.TH1D("reslY","Laser track residual Y", 200,-maxDist, maxDist)
resY.GetXaxis().SetTitle("cm")
reslZ = ROOT.TH1D("reslZ","Laser track residual Z", 200,-maxDist, maxDist)
resZ.GetXaxis().SetTitle("cm")

dr3D = ROOT.TH3D("dr3D","Radial residuals",500,-39.5,109.5,200,15.5,41.5,300,-2,2)




rfile = ROOT.TFile(tfile, "read")
tree = rfile.Get("cbmsim")
tree.SetBranchStatus("*", 0)
tree.SetBranchStatus("TpcLaserMatchingResidual.*", 1)
tree.SetBranchStatus("TpcGFLaserTrack.*", 1)
tree.SetBranchStatus("TpcTrackResiduals.*", 1)


for e in tree:
    for res in e.TpcLaserMatchingResidual :
        residual = res.getResXYZ()
        resX.Fill(residual.X())
        resY.Fill(residual.Y())
        resZ.Fill(residual.Z())
    count = 0    
    for rc in e.TpcTrackResiduals :
        trackID = rc.getRefIndex()
        #get corresponding track:
        track = e.TpcGFLaserTrack.At(trackID)
        for ri in range(rc.getSize()) :
            #track Direction:
            refPlane = track.getCardinalRep().getReferencePlane()  #workaround
            tDir = track.getCardinalRep().getMom(refPlane)
            tDir.SetMag(1.)
            res3 = rc.getResidual(ri).getResXYZ()
            reslX.Fill(res3.X())
            reslY.Fill(res3.Y())
            reslZ.Fill(res3.Z())
            pP = rc.getResidual(ri).getRefPos()
            pDir = ROOT.TVector3(pP)   #copy
            pDir.SetZ(0.)
            pDir.SetMag(1.)
            dr = pDir*res3
            dr3D.Fill(pP.Z(),pP.Perp(),dr)
            
            
canv = ROOT.TCanvas("blah", "Matching Residuals", 400,1000)
canv.Divide(1,3)
canv.cd(1)
resX.Draw()
canv.cd(2)
resY.Draw()
canv.cd(3)
resZ.Draw()

canv2 = ROOT.TCanvas("blah2", "Track Residuals", 400,1000)
canv2.Divide(1,3)
canv2.cd(1)
reslX.Draw()
canv2.cd(2)
reslY.Draw()
canv2.cd(3)
reslZ.Draw()

canv3 = ROOT.TCanvas()
dr3D.Draw()


       
    
input()