#do python /scratch/sdorheim/FOPI/newDev/OUTPUT/git_scripts//Alignment/runStandaloneResCalc.py -n 0 -r $i -m 3d -alFile alinput.txt -p /nfs/mds/data/tpc/fopi/2011/reconstructed/smooth_finalCandidate3_newAlignment_newZshift/ -fb TpcFopiRecoSGE -smoothing -field -pt 0.15 -max-pt 1.5 -sectorCut 0 -fPath /nfs/mds/data/tpc/fopi/2011/fopi_reco_sverre/rev176/; done


import ROOT, glob, math, sys, os, argparse
print "####################################################################################################################################################################################"
print sys.argv[0] + "running in this directory:"
os.system("pwd")
print "####################################################################################################################################################################################"
parser=argparse.ArgumentParser(description="Scrip to run one Alignment step/iteraton")
parser.add_argument("-n","--number-of-events", help="number of input events to be used for alignment", type=int, default=140000)
parser.add_argument("-r","--run", help="run number to use if n<140k otherwhise specify run range 32XX , 34XX, 35XX", type=int, default=3280)
parser.add_argument("-nr","--n-runs",help="number of runs", type=int, default=1)
parser.add_argument("-m", "--method", help="residual calculator to be used, \n 2D, 3D, 3D-simple, track-track, cdc", type=str, default="3d")
parser.add_argument("-d", "--dummy", help="enable dummy mode, print options but do not loop", action="store_true")
parser.add_argument("-alFile", "--alignmentfile", help="specify alignment file to be used", type=str, default="alinput.txt")
parser.add_argument("-p", "--path", help="Specify path to reconstructed files", type=str, default = "/scratch/other/sdorheim/DATA/unSmooth_noV0FopiTracking_RpcBarPixFits_DAF/")
parser.add_argument("-fb", "--file-base", help="Specify file base of reco files. Used as follows: \n \"path + filebase+\"_run\"+runStr+\"_part\"+partStr+\"*\"+smoothingText+\".reco.root\"", type=str, default="TpcFopiRecoPersist_run")
parser.add_argument("-smoothing","--smoothing", help="Flag if using smoothed data", action="store_true")
parser.add_argument("-charge","--charge-filter", help="Specify filtering on charge, =0: no filtering, >0: only positive charge, <0: only negative charge", type=int, default=0)
parser.add_argument("-field", "--mag-field", help="flag to turn on magnetic field", action="store_true")
parser.add_argument("-pt", "--pt-limit", help="limit on Cdc track pt>limit", default=0.8, type=float)
parser.add_argument("-max-pt", "--maximum-pt", help="Maximum Cdc track pt, pt<limit", default=99, type=float)
parser.add_argument("-sectorCut", "--sectorCut", help="Cut on n sectors", default=0, type=int)
parser.add_argument("-cdcBranch", "--cdc-branch-name", help="cdc branch to be used", default="CdcTrackPostFit")
parser.add_argument("-t0suffix", "--t0-suffix", help="t0 suffix for t0-scan, exampel: _t0_-197.61399841",default="")
parser.add_argument("-fPath","--fopi-path", help="location of fopi-reco files", default = "/nfs/mds/data/tpc/fopi/2011/fopi_reco_sverre/rev171")
parser.add_argument("-cdc-mass","--cdc-mass-cut", help="only consider tracks with mass smaller thatn..", default = -1., type=float)
parser.add_argument("-outfileprefix", "--outFilePrefix", help="name of output files will be <prefix>+_+<method>+_<run>+.root", default="residuals")
parser.add_argument("-storeHistos", "--storeHistos", help="store histos", action="store_true")
parser.add_argument("-cutOnRpc", "--cutOnRpc", help="cut out tracks without rpc hit", action="store_true")
parser.add_argument("-minCdcHits", "--minCdcHits", help=" minimum number of cdc hits needed", type=int, default=45)
args=parser.parse_args()
print sys.argv[0], "called with these args"


for arg, value in vars(args).items():
    print arg,": ", value

dummy=args.dummy
#argument parsing:

if args.smoothing:
    smoothingText="smoothed"
else:
    smoothingText=""

method=args.method

if not (method.lower() == "2d" or method.lower() == "3d" or method.lower() == "3d-simple" or method.lower() == "track-track" or method.lower() == "cdc" or method.lower()=="tpcmc"):
    print "Allowed methods: 2D, 3D, 3D-simple, track-track, cdc, tpcmc"
    print "you used: ", method.lower()
    exit()


ROOT.gROOT.ProcessLine(".x $VMCWORKDIR/gconfig/rootlogon.C")
os.system("echo $VMCWORKDIR")

timer=ROOT.TStopwatch()
timer.Start(True)

runStr=str(args.run)

TpcSAResCalc=ROOT.TpcStandaloneResCalc()
TpcSAResCalc.setAlignmentFile(args.alignmentfile)
TpcSAResCalc.setOutFileName(args.outFilePrefix+"_"+method.lower()+"_"+str(args.run)+".root")
TpcSAResCalc.setHistoFileName("AlignmentHist_"+method.lower()+"_"+str(args.run)+".root")

partStr=""

#needed for parameter managment
fRun=ROOT.FairRunAna()
ROOT.gROOT.SetBatch()

chain=ROOT.TChain("cbmsim")
t0suffix=args.t0_suffix
if t0suffix=="\'\'":
    t0suffix=""
print t0suffix

if args.run !="0":
    print "test"
    runStart=args.run
    files=[]
    

    for run in range(runStart,runStart+1):
        searchstring= args.path + args.file_base+str(run)+t0suffix+"_part*"+smoothingText+".reco.root"
        runFiles=glob.glob(searchstring)
        files.extend(runFiles)

    logfile=open("runStandaloneResCalc_"+method.lower()+"_"+str(args.run)+".txt",'w')
    logfile.write("Alignment run with these parameters\n")
    for arg, value in vars(args).items():
        logfile.write( str(arg)+": "+ str(value)+"\n")
    logfile.write("files where searched using this searchstring:\n "+searchstring)
    logfile.close()
    print "searching for files fitting to this search string:\n"+searchstring
    
    files.sort()
    
    print len(files), " files found"
    
    firstRun=True
    for fileN in files :
        print fileN
        #    lambdaChain.addFile(file)
        if firstRun:
            firstFileN=fileN
            firstRun=False
        chain.AddFile(fileN)

searchString=args.fopi_path+"/*"+str(args.run)+"*.root"
fopiFiles=glob.glob(searchString)
print args.t0_suffix
print fopiFiles
if len(fopiFiles)!= 1:
    print "more than one fopi reco file found, sort it out:"
    print searchString
    print fopiFiles
else:
    fopiFile=ROOT.TFile(fopiFiles[0],"read")
    print fopiFiles
    fopiChain=fopiFile.Get("FOPITree")
    TpcSAResCalc.setFopiChain(fopiChain)
    
fRun.SetInputFile(firstFileN)
#chain.SetBranchStatus("*.*",0);

chain.SetCacheLearnEntries(1)
chain.SetCacheSize(1000000)
print "\n"
print "chain entries:"
print chain.GetEntries()
print "\n"

if chain.GetEntries<1:
    exit()
if args.mag_field:
    print "\n\n\n\nUSING FOPI field\n\n\n\n"
    fMagField = ROOT.FOPIField(0.616); 
    fMagField.SetTargetOffset(-61);
else:
    print "\n\n\n\nUSING const field\n\n\n\n"
    fMagField=ROOT.PndConstField();
    fMagField.SetField(0., 0. , 0.616 ); #values are in kG
    fMagField.SetFieldRegion(-120,120,-120, 120, -2000, 2000); # values are in cm
fRun.SetField(fMagField)

ROOT.GFFieldManager.getInstance().init(ROOT.PndFieldAdaptor(fRun.GetField()));
ROOT.GFMaterialEffects.getInstance().init(ROOT.GFTGeoMaterialInterface());
################################################################################

# check is the geometry loaded via the param file
taroff=-61 
targetpos=taroff+41
targetmat="carbon"
#geoFile = ROOT.TString()
geoFile="/tmp/FopiGeom_s339"+str(targetpos)+str(args.run)+".root"
#print cmd
cmd=" root -b -q '$FOPI2ROOT/fopigeometry/FopiGeom.C("+str(targetpos)+",\""+targetmat+"\",\""+args.alignmentfile+"\",\""+geoFile+"\",kFALSE)'"
os.system(cmd)
print "################################################################################"
print geoFile
print "################################################################################"

fRun.SetGeomFile(geoFile);

#does this work on real data?
basedir = ROOT.gSystem.Getenv("VMCWORKDIR");
rtdb = fRun.GetRuntimeDb();
parInput1=ROOT.FairParRootFileIo()

parInput1.open(firstFileN.replace("reco.root","param.root"))
rtdb.setFirstInput(parInput1);

fRun.SetOutputFile("dummy"+runStr+".root");


nEv=args.number_of_events
if nEv==0:
    nEv=0

print "\n\n\n\n\n\n\n\n\n#### fRun.Init() ################################################################################\n\n"
tpcDrifter = ROOT.TpcDriftTask()
tpcDrifter.SetDistort(False)
fRun.AddTask(tpcDrifter)
fRun.Init()
print "\n\n\n\n\n\n\n\n\n#### fRun.Init() done ###########################################################################\n\n"

TpcSAResCalc.setMaxEvents(nEv)
print "TpcSAResCalc.setMaxEvents(nEv)"
#set limits on histograms

print "TpcSAResCalc.setInputChain(chain)"
TpcSAResCalc.setInputChain(chain)
print "set det name"
TpcSAResCalc.setDetectorName("tpc")
resCalc=None
if method.lower() == "3d":
    resCalc=ROOT.CdcTrackTpcHitResCalc()
    resCalc.setHistoFileName("CdcTrackTpcHitResCalcQA"+"_"+str(args.run)+".root")
    #    resCalc.addBranchName("TpcCluster")
    resCalc.addBranchName(args.cdc_branch_name)
    resCalc.addBranchName("TpcTrackPostFit")
    resCalc.addBranchName("TpcSPHit")
    resCalc.addBranchName("FopiTrackTuples")
    #   resCalc.addBranchName("FopiTrackPostFit")
    resCalc.setUseFelixMatching("FopiTrackTuples")
    
    TpcSAResCalc.addBranch(args.cdc_branch_name,"GFTrack")
    TpcSAResCalc.addBranch("TpcTrackPostFit","GFTrack")
    TpcSAResCalc.addBranch("TpcSPHit","TpcSPHit")
    TpcSAResCalc.addBranch("FopiTrackTuples", "MatchingTuple")
    #    TpcSAResCalc.addBranch("TpcCdcMatchingPair", "MatchingTuple")
    #    TpcSAResCalc.addBranch("FopiTrackPostFit", "GFTrack")
    
    resCalc.setCheckFopiTrack(False)
    #    TpcSAResCalc.addBranch("FopiTrackTuples", "MatchingTuple")
    #    TpcSAResCalc.addBranch("TpcCdcMatchingPair", "MatchingPair")
    #    TpcSAResCalc.addBranch("FopiTrackPostFit", "GFTrack")

    resCalc.setTpcClusterBranchName("TpcCluster")
    resCalc.setCdcGFTrackBranchName(args.cdc_branch_name)
    resCalc.setTpcTrackBranchName("TpcTrackPostFit")
    resCalc.setTpcSPHitBranchName("TpcSPHit")
    
    resCalc.setTpcPtCut(0.0) # might defaulted to something other than 0, since tpc has ~no pt-resolution alone, this cut makes no sense
    resCalc.setCdcPtCut(args.pt_limit)
    resCalc.setCdcMaxPtCut(args.maximum_pt)
    resCalc.setThetaCut(20)
    resCalc.setCdcExtrapolateToZ(True)
    resCalc.setCutOnRpc(args.cutOnRpc)    
    resCalc.setThetaCut(20)
    #    resCalc.setUseClusters(True)
    resCalc.setOrigTransfName("orig")
    if args.cdc_mass_cut>0:
        resCalc.setMassCut(args.cdc_mass_cut)
    resCalc.setMinHitsCdc(args.minCdcHits)
    #    resCalc.applyDevMap("/scratch/sdorheim/FOPI/newDev/tpc/FOPI/par/full_field_-148.62_-236.74_4s3_4.700_0.000sh0.000_65ar81_merged.txt")
    #    resCalc.applyDevMap("/scratch/sdorheim/FOPI/newDev/tpc/FOPI/par/full_field_-148.62_-232.79_4s3_4.700_-0.000sh0.000_0.0000_65ar80__merged.txt")
    if args.sectorCut>0 or args.cdc_mass_cut>0:
        
        resCalc.addBranchName("CdcHitBr")
        resCalc.addBranchName("CdcTrackBr")
        resCalc.setCdcTrackBranchName("CdcTrackBr")
        resCalc.setCdcHitBranchName("CdcHitBr")
        resCalc.setCutOnNSectors(args.sectorCut)
        TpcSAResCalc.addFopiBranch("CdcHitBr","CdcHit")
        TpcSAResCalc.addFopiBranch("CdcTrackBr","CdcTrack")
        TpcSAResCalc.setTpcEventBranchName("TpcEventIdentifier")
        TpcSAResCalc.setFopiEventBranchName("FopiEventBr")
    #        resCalc.setUseClusters(True)
    residualcCalculatorName="CdcTrackTpcHitResCalc"

    if True:    
        resCalc.setCheckFopiTrack(False)
if args.storeHistos:
    TpcSAResCalc.setSaveHistos(True)       
TpcSAResCalc.setResCalculatior(resCalc,residualcCalculatorName)
TpcSAResCalc.setChargeFiltering(args.charge_filter)

TpcSAResCalc.setMaxRangeXY(2);

TpcSAResCalc.init()
print "###################################################################################################"
print "###################################################################################################"
print "###################################################################################################"
print "###################################################################################################"
print "# "
print "# "
#print  "# ",ROOT.TpcDigiMapper.getInstance().getGas().VDrift(0,0);
print "# "
print "# "
print "###################################################################################################"
print "###################################################################################################"
print "###################################################################################################"
print "###################################################################################################"

timer.Stop();
rtime = timer.RealTime()
ctime = timer.CpuTime()

timer.Start(True)
print "residual preparation took ", rtime, "and ctime: ",ctime
TpcSAResCalc.execute()
timer.Stop()
rtime2 = timer.RealTime()
ctime2 = timer.CpuTime()
print "residual preparation took ", rtime, "and ctime: ",ctime
print "residual calculation took ", rtime2, "and ctime: ",ctime2
if args.storeHistos:
    TpcSAResCalc.writeHistos()

print "\n\n\n\n done \n\n\n\n "
if args.storeHistos:
    if method.lower()=="3d" or method.lower()=="2d" or method.lower() == "cdc" and not mc:
        resCalc.writeHists()
print "\n\n\n\n done \n\n\n\n "
exit()