Python RooArgList.getSize方法代码示例

# 需要导入模块: from ROOT import RooArgList [as 别名]
# 或者: from ROOT.RooArgList import getSize [as 别名]
def findOnePe(hist, ws, name='x', Npe = 1):
    fitPed(hist, ws, name)
    x = ws.var(name)

    ped = ws.pdf('ped')
    pedWidth = ws.var('pedWidth')

    pdfs = RooArgList(ped)
    pdfList = []

    fped = RooRealVar('fped', 'f_{ped}', 0.8, 0., 1.)
    fractions = RooArgList(fped)
    fList = []
    peList = []

    peMean = RooRealVar('peMean', 'mean_{pe}', 6., 0., 20.)
    peWidth = RooRealVar('peWidth', 'width_{pe}', pedWidth.getVal(), 0., 10.)

    for i in range(0, Npe):
        pem = RooFormulaVar('pem{0}'.format(i+1), '@0+{0}*@1'.format(i+1),
                            RooArgList(ws.var('pedMean'), peMean))
        peList.append(pem)
        npepdf = RooGaussian('pe{0}pdf'.format(i+1), 'pe{0}pdf'.format(i+1),
                             x, pem, pedWidth)
        pdfs.add(npepdf)
        pdfList.append(npepdf)

        fnpe = RooRealVar('fpe{0}'.format(i+1), 'fpe{0}'.format(i+1),
                          0.5, -0.1, 1.0)
        fractions.add(fnpe)
        fList.append(fnpe)

    #bgMean = RooRealVar("bgMean", "bgMean", 6.0, x.getMin(), x.getMax())
    bgScale = RooRealVar("bgScale", "bgScale", 0.5, -1.0, Npe + 1.0)
    bgMean = RooFormulaVar("bgMean", "@[email protected]*@2",
                           RooArgList(peMean, ws.var('pedMean'), bgScale))
    bgWidthL = RooRealVar("bgWidthL", "bgWidthL", pedWidth.getVal()*2,
                          0., 25.)
    bgWidthR = RooRealVar("bgWidthR", "bgWidthR", pedWidth.getVal()*7,
                          0., 25.)

    bgGauss = RooBifurGauss("bgGauss", "bgGauss", x, bgMean,
                            bgWidthR, bgWidthR)

    if (Npe > 1):
        pdfs.add(bgGauss)
    else:
        fractions.remove(fractions.at(fractions.getSize()-1))

##     pem = RooFormulaVar('pem', '@[email protected]', RooArgList(peMean, ws.var('pedMean')))
##     firstPe = RooGaussian('firstPe', 'firstPe', x, pem, peWidth)

##     pdfs.Print("v")
##     fractions.Print("v")
    pedPlusOne = RooAddPdf('pedPlusOne', 'pedPlusOne', pdfs, fractions, True)

    ## pedWidth = ped.GetParameter(2)
    ## pedMean = ped.GetParameter(1)
    ## pedA = ped.GetParameter(0)
    
    secondMax = hist.GetMaximumBin() + 1
    goingDown = True
    maxVal = hist.GetBinContent(secondMax)
    foundMax = False
    while (not foundMax) and (secondMax < hist.GetNbinsX()):
        tmpVal = hist.GetBinContent(secondMax+1)
        if (tmpVal < maxVal):
            if not goingDown:
                foundMax = True
            else:
                goingDown = True
                maxVal = tmpVal
                secondMax += 1
        elif (tmpVal > maxVal):
            goingDown = False
            maxVal = tmpVal
            secondMax += 1
        else:
            maxVal = tmpVal
            secondMax += 1

    secondMaxx = hist.GetBinCenter(secondMax)
    print 'found 2nd maximum in bin',secondMax,'value',secondMaxx

##     peMean.setVal(secondMaxx)
##     bgMean.setVal(secondMaxx*0.6)
    x.setRange('pedPlus_fit', x.getMin(), ws.var('pedMean').getVal()+pedWidth.getVal()*6.*(Npe+0))

    pedPlusOne.fitTo(ws.data('ds'), RooFit.Minos(False),
                     RooFit.Range('pedPlus_fit'),
                     RooFit.PrintLevel(1))

    getattr(ws, 'import')(pedPlusOne)

# 需要导入模块: from ROOT import RooArgList [as 别名]
# 或者: from ROOT.RooArgList import getSize [as 别名]
# globals, and aliases
from ROOT import gPad, gStyle, gSystem, kRed, kBlue, kAzure, kGreen, kBlack
const = RooRealConstant.value

# setup
gStyle.SetOptStat('nemrou')
canvas = TCanvas('canvas', 'canvas', 800, 600)
canvas.Print('{}['.format(plotfile))

# read workspace
ffile = TFile.Open(rfile, 'read')
workspace = ffile.Get('workspace')

pdfs = RooArgList(workspace.allPdfs())
for i in range(pdfs.getSize()):
    name = pdfs[i].GetName()
    if name.find(mode) >= 0:
        mykpdf = pdfs[i]
assert(mykpdf)

## variables
time = RooRealVar('time', 'Time [ps]', 0.2, 15.0)
time.setBins(bins)
time.setBins(bins*3, 'cache')

kfactor = workspace.var('kfactorVar')
kfactor.setRange(0.85, 1.05)
gamma = RooRealVar('gamma', 'gamma', 0.661, 0., 3.)
kgamma = RooProduct('kgamma', 'kgamma', RooArgList(gamma, kfactor))
dGamma = RooRealVar('dGamma', 'dGamma', -0.106, -3., 3.)

# 需要导入模块: from ROOT import RooArgList [as 别名]
# 或者: from ROOT.RooArgList import getSize [as 别名]
    parlist = fitresult.floatParsFinal()
    cmatrix = fitresult.covarianceMatrix()

    veclist = RooArgList()
    for i in range(parlist.getSize()):
        name = "%s_%d" % (parlist[i].GetName(), i)
        veclist.add(parlist[i].clone(name))

    multigauss = RooMultiVarGaussian("multigauss", "multigauss", veclist, parlist, cmatrix)
    dset = multigauss.generate(RooArgSet(veclist), 1000)

    fns = []
    for entry in range(dset.numEntries()):
        vecset = dset.get(entry)
        veclist = RooArgList(vecset)
        for pars in range(veclist.getSize()):
            acceptance.SetParameter(pars, veclist[pars].getVal())
        fns += [acceptance.Clone("%s_%d" % (acceptance.GetName(), entry))]

    avgfn = BinnedAvgFunction(fns, xbincs)
    avgfn.calculate()
    accfns += [avgfn.get_avg_fn()]
    accfnerrs += [avgfn.get_avg_fn_var()]


means = numpy.zeros(nbins, dtype=float)
varis = numpy.zeros(nbins, dtype=float)
for ibin in range(nbins):
    means[ibin] = accfns[0][ibin] / accfns[1][ibin]
    varis[ibin] = accfnerrs[0][ibin] + accfnerrs[1][ibin]

# 需要导入模块: from ROOT import RooArgList [as 别名]
# 或者: from ROOT.RooArgList import getSize [as 别名]
usig2 = 0.
totalYield = 0.
sigYield = 0.
sigErrs = {}

sigYieldFilename = 'last_H%i_%s_%iJets_signalYield.txt' % (opts.mH, 
                                                           modeString, opts.Nj)
sigYieldsFile = open(sigYieldFilename, 'w')

WpJNonPoissonError = 0

print
print '-------------------------------'
print 'Yields in signal box'
print '-------------------------------'
for i in range(0, yields.getSize()):
    theName = yields.at(i).GetName()
    if theName[0] == 'n':
        totalYield += yields.at(i).getVal()
        theIntegral = 1.
        if (theName == 'nDiboson'):
            theIntegral = dibosonInt.getVal()/dibosonFullInt.getVal()
        elif (theName == 'nWjets'):
            theIntegral = WpJInt.getVal()/WpJFullInt.getVal()
            if (yields.at(i).getError()**2 > yields.at(i).getVal()):
                WpJNonPoissonError = sqrt(yields.at(i).getError()**2 - \
                                              yields.at(i).getVal())
            else:
                WpJNonPoissonError = 0.
        elif (theName == 'nTTbar'):
            theIntegral = ttbarInt.getVal()/ttbarFullInt.getVal()

# 需要导入模块: from ROOT import RooArgList [as 别名]
# 或者: from ROOT.RooArgList import getSize [as 别名]
def accbuilder(time, knots, coeffs):
    # build acceptance function
    from copy import deepcopy
    myknots = deepcopy(knots)
    mycoeffs = deepcopy(coeffs)
    from ROOT import (RooBinning, RooArgList, RooPolyVar,
            RooCubicSplineFun)
    if (len(myknots) != len(mycoeffs) or 0 >= min(len(myknots), len(mycoeffs))):
        raise ValueError('ERROR: Spline knot position list and/or coefficient'
                'list mismatch')
    # create the knot binning
    knotbinning = WS(ws, RooBinning(time.getMin(), time.getMax(), 'knotbinning'))
    for v in myknots:
        knotbinning.addBoundary(v)
    knotbinning.removeBoundary(time.getMin())
    knotbinning.removeBoundary(time.getMax())
    knotbinning.removeBoundary(time.getMin())
    knotbinning.removeBoundary(time.getMax())
    oldbinning, lo, hi = time.getBinning(), time.getMin(), time.getMax()
    time.setBinning(knotbinning, 'knotbinning')
    time.setBinning(oldbinning)
    time.setRange(lo, hi)
    del knotbinning
    del oldbinning
    del lo
    del hi
    # create the knot coefficients
    coefflist = RooArgList()
    i = 0
    for v in mycoeffs:
        coefflist.add(WS(ws, RooRealVar('SplineAccCoeff%u' % i,
            'SplineAccCoeff%u' % i, v)))
        i = i + 1
    del mycoeffs
    coefflist.add(one)
    i = i + 1
    myknots.append(time.getMax())
    myknots.reverse()
    fudge = (myknots[0] - myknots[1]) / (myknots[2] - myknots[1])
    lastmycoeffs = RooArgList(
            WS(ws, RooConstVar('SplineAccCoeff%u_coeff0' % i,
                'SplineAccCoeff%u_coeff0' % i, 1. - fudge)),
            WS(ws, RooConstVar('SplineAccCoeff%u_coeff1' % i,
                'SplineAccCoeff%u_coeff1' % i, fudge)))
    del myknots
    coefflist.add(WS(ws, RooPolyVar(
        'SplineAccCoeff%u' % i, 'SplineAccCoeff%u' % i,
        coefflist.at(coefflist.getSize() - 2), lastmycoeffs)))
    del i
    # create the spline itself
    tacc = WS(ws, RooCubicSplineFun('SplineAcceptance', 'SplineAcceptance', time,
        'knotbinning', coefflist))
    del lastmycoeffs
    # make sure the acceptance is <= 1 for generation
    m = max([coefflist.at(j).getVal() for j in
        xrange(0, coefflist.getSize())])
    from ROOT import RooProduct
    c = WS(ws, RooConstVar('SplineAccNormCoeff', 'SplineAccNormCoeff', 0.99 / m))
    tacc_norm = WS(ws, RooProduct('SplineAcceptanceNormalised',
        'SplineAcceptanceNormalised', RooArgList(tacc, c)))
    del c
    del m
    del coefflist
    return tacc, tacc_norm

# 需要导入模块: from ROOT import RooArgList [as 别名]
# 或者: from ROOT.RooArgList import getSize [as 别名]
from ROOT import gPad, gStyle, kRed, kBlue, kAzure, kGreen, kBlack

gStyle.SetOptStat('nemrou')
canvas = TCanvas('canvas', 'canvas', 800, 600)
canvas.Print('{}['.format(plotfile))

ffile = TFile.Open(rfile, 'read')

if ifpdf:
    workspace = ffile.Get('workspace')

    kfactor = workspace.var('kfactorVar')
    kfactor.setRange(0.9, 1.1)
    pdfs = RooArgList(workspace.allPdfs())

    for i in range(pdfs.getSize()):
        fr = kfactor.frame()
        pdfs[i].plotOn(fr, RooFit.LineColor(kBlack), RooFit.FillColor(kAzure+1),
                       RooFit.DrawOption('lf')) # FIXME: doesn't draw the line!'
        fr.Draw()
        canvas.Print(plotfile)
else:
    modes = {}
    klist = ffile.GetListOfKeys()
    for item in klist:
        name = item.GetName()
        if not name.startswith('mBresn'):
            print 'MSG: Skipping, unknown object: %s' % name
            continue # ntuples are named mBresn_*
        sample = name.split('_')
        modes[sample[1]] = modes.get(sample[1],[]) + [item] # up and down for each mode

# 需要导入模块: from ROOT import RooArgList [as 别名]
# 或者: from ROOT.RooArgList import getSize [as 别名]

#.........这里部分代码省略.........
                'DeltaGamma': DeltaGamma,
                'DeltaM': DeltaM,
                'C': C, 'D': D, 'Dbar':Dbar, 'S': S, 'Sbar': Sbar,
                'timeresmodel': timeresmodel,
                'acceptance': acceptance,
                'timeerrpdf': timeerrpdf,
                'mistagpdf': mistagpdf,
                'mistagobs': mistagobs,
                'kfactorpdf': kfactorpdf,
                'kvar': kvar,
                'aprod': aprod,
                'adet': adet,
                'atageff': atageff
                }
        print 'buildBDecayTimePdf('
        for kw in kwargs:
            print '\t%s = %s' % (kw, kwargs[kw])
        print '\t)'
        print 72 * '#'

    # constants used
    zero = WS(ws, RooConstVar('zero', 'zero', 0.))
    one = WS(ws, RooConstVar('one', 'one', 1.))

    if None == aprod: aprod = zero
    if None == adet: adet = zero
    if None == atageff: atageff = [ zero ]
    if None == mistagpdf:
        mistagobs = None
    else: # None != mistagpdf
        if None == mistagobs:
            raise NameError('mistag pdf set, but no mistag observable given')

    # if no time resolution model is set, fake one
    if timeresmodel == None:
        timeresmodel = WS(ws, RooTruthModel('%s_TimeResModel' % name,
            '%s time resolution model' % name, time))
    elif timeresmodel == 'Gaussian':
        timeresmodel = WS(ws, RooGaussModel('%s_TimeResModel' % name,
            '%s time resolution model' % name, time, zero, timeerr))

    # apply acceptance (if needed)
    timeresmodel = applyBinnedAcceptance(
            config, ws, time, timeresmodel, acceptance)
    if config['UseKFactor']:
        timeresmodel = applyKFactorSmearing(config, ws, time, timeresmodel,
                kvar, kfactorpdf, [ Gamma, DeltaGamma, DeltaM ])
    if config['ParameteriseIntegral']:
        parameteriseResModelIntegrals(config, ws, timeerrpdf, timeerr, timeresmodel)

    # if there is a per-event mistag distributions and we need to do things
    # correctly
    if None != mistagpdf:
        otherargs = [ mistagobs, RooArgList(*mistagpdf), RooArgList(*tageff) ]
    else:
        otherargs = [ RooArgList(*tageff) ]
    bcalib = RooArgList()
    bbarcalib = RooArgList()
    for t in mistag:
        bcalib.add(t[0])
        if len(t) > 1:
            bbarcalib.add(t[1])
    otherargs.append(bcalib)
    if (bbarcalib.getSize()):
        otherargs.append(bbarcalib)
    otherargs += [ aprod, adet, RooArgList(*atageff) ]
    flag = 0
    if 'Bs2DsK' == name and 'CADDADS' == config['Bs2DsKCPObs']:
        flag = DecRateCoeff.AvgDelta
    # build coefficients to go into RooBDecay
    cosh = WS(ws, DecRateCoeff('%s_cosh' % name, '%s_cosh' % name,
        DecRateCoeff.CPEven, qf, qt, one, one, *otherargs))
    sinh = WS(ws, DecRateCoeff('%s_sinh' % name, '%s_sinh' % name,
        flag | DecRateCoeff.CPEven, qf, qt, D, Dbar, *otherargs))
    cos = WS(ws, DecRateCoeff('%s_cos' % name, '%s_cos' % name,
        DecRateCoeff.CPOdd, qf, qt, C, C, *otherargs))
    sin = WS(ws, DecRateCoeff('%s_sin' % name, '%s_sin' % name,
        flag | DecRateCoeff.CPOdd | DecRateCoeff.Minus,
        qf, qt, S, Sbar, *otherargs))
    del flag
    del otherargs

    # build (raw) time pdf
    tau = WS(ws, Inverse('%sTau' % Gamma.GetName(),
        '%s #tau' % Gamma.GetName(), Gamma))
    retVal = WS(ws, RooBDecay(
        '%s_RawTimePdf' % name, '%s raw time pdf' % name,
        time, tau, DeltaGamma,  cosh, sinh, cos, sin,
        DeltaM, timeresmodel, RooBDecay.SingleSided))

    retVal = applyDecayTimeErrPdf(config, name, ws, time, timeerr, qt, qf,
            mistagobs, retVal, timeerrpdf, mistagpdf)
    
    # if we do not bin the acceptance, we apply it here
    retVal = applyUnbinnedAcceptance(config, name, ws, retVal, acceptance)

    retVal.SetNameTitle('%s_TimePdf' % name, '%s full time pdf' % name)

    # return the copy of retVal which is inside the workspace
    return WS(ws, retVal)

# 需要导入模块: from ROOT import RooArgList [as 别名]
# 或者: from ROOT.RooArgList import getSize [as 别名]

#.........这里部分代码省略.........
                                args.add(iset.find(n))
                    var = WS(fws, RooFormulaVar(name, name, name, args))
                    addiset.addClone(var)
            else:
                for dsname in ((config['DataSetNames'], )
                        if type(config['DataSetNames']) == str else
                        config['DataSetNames']):
                    break
                leaf = f.Get(dsname).GetLeaf(name)
                if None == leaf:
                    leaf = f.Get(dsname).GetLeaf(name + '_idx')
                if leaf.GetTypeName() in (
                        'char', 'unsigned char', 'Char_t', 'UChar_t',
                        'short', 'unsigned short', 'Short_t', 'UShort_t',
                        'int', 'unsigned', 'unsigned int', 'Int_t', 'UInt_t',
                        'long', 'unsigned long', 'Long_t', 'ULong_t',
                        'Long64_t', 'ULong64_t', 'long long',
                        'unsigned long long'):
                    var = WS(fws, RooCategory(name, name))
                    tit = obj.typeIterator()
                    ROOT.SetOwnership(tit, True)
                    while True:
                        tobj = tit.Next()
                        if None == tobj: break
                        var.defineType(tobj.GetName(), tobj.getVal())
                else:
                    var = WS(fws, RooRealVar(name, name, -sys.float_info.max,
                        sys.float_info.max))
                iset.addClone(var)
        for dsname in ((config['DataSetNames'], )
               if type(config['DataSetNames']) == str else
               config['DataSetNames']):
            tmpds = WS(fws, RooDataSet(dsname, dsname,f.Get(dsname), iset), [])
            if 0 != addiset.getSize():
                # need to add columns with RooFormulaVars
                tmpds.addColumns(addiset)
            del tmpds
    # local data conversion routine
    def doIt(config, rangeName, dsname, sname, names, dmap, dset, ddata, fws):
        sdata = fws.obj(dsname)
        if None == sdata: return 0
        if None != config['DataSetCuts']:
            # apply any user-supplied cuts
            newsdata = sdata.reduce(config['DataSetCuts'])
            ROOT.SetOwnership(newsdata, True)
            del sdata
            sdata = newsdata
            del newsdata
        sset = sdata.get()
        smap = { }
        for k in names:
            smap[k] = sset.find(config['DataSetVarNameMapping'][k])
        if 'sample' in smap.keys() and None == smap['sample'] and None != sname:
            smap.pop('sample')
            dmap['sample'].setLabel(sname)
        if None in smap.values():
            raise NameError('Some variables not found in source: %s' % str(smap))
        # # additional complication: toys save decay time in ps, data is in nm
        # # figure out which time conversion factor to use
        # timeConvFactor = 1e9 / 2.99792458e8
        # meantime = sdata.mean(smap['time'])
        # if ((dmap['time'].getMin() <= meantime and
        #         meantime <= dmap['time'].getMax() and config['IsToy']) or
        #         not config['IsToy']):
        #     timeConvFactor = 1.
        # print 'DEBUG: Importing data sample meantime = %f, timeConvFactor = %f' % (

# 需要导入模块: from ROOT import RooArgList [as 别名]
# 或者: from ROOT.RooArgList import getSize [as 别名]
def buildSplineAcceptance(
        ws,     # workspace into which to import
        time,   # time variable
        pfx,    # prefix to be used in names
        knots,  # knots
        coeffs, # acceptance coefficients
        floatParams = False, # float acceptance parameters
        debug = False # debug printout
        ): 
    """
    build a spline acceptance function

    ws          -- workspace into which to import acceptance functions
    time        -- time observable
    pfx         -- prefix (mode name) from which to build object names
    knots       -- list of knot positions
    coeffs      -- spline coefficients
    floatParams -- if True, spline acceptance parameters will be floated
    debug       -- if True, print some debugging output

    returns a pair of acceptance functions, first the unnormalised one for
    fitting, then the normalised one for generation

    The minimum and maximum of the range of the time variable implicitly
    defines the position of the first and last knot. The other knot positions
    are passed in knots. Conversely, the coeffs parameter records the height
    of the sline at all but the last two knot positions. The next to last knot
    coefficient is fixed to 1.0, thus fixing the overall scale of the
    acceptance function. The spline coefficient for the last knot is fixed by
    extrapolating linearly from the two knots before; this prevents
    statistical fluctuations at the low stats high lifetime end of the
    spectrum to curve the spline.
    """
    # build acceptance function
    from copy import deepcopy
    myknots = deepcopy(knots)
    mycoeffs = deepcopy(coeffs)
    from ROOT import (RooBinning, RooArgList, RooPolyVar, RooCubicSplineFun,
            RooConstVar, RooProduct, RooRealVar)
    if (len(myknots) != len(mycoeffs) or 0 >= min(len(myknots), len(mycoeffs))):
        raise ValueError('ERROR: Spline knot position list and/or coefficient'
                'list mismatch')
    one = WS(ws, RooConstVar('one', '1', 1.0))
    # create the knot binning
    knotbinning = WS(ws, RooBinning(time.getMin(), time.getMax(),
        '%s_knotbinning' % pfx))
    for v in myknots:
        knotbinning.addBoundary(v)
    knotbinning.removeBoundary(time.getMin())
    knotbinning.removeBoundary(time.getMax())
    knotbinning.removeBoundary(time.getMin())
    knotbinning.removeBoundary(time.getMax())
    oldbinning, lo, hi = time.getBinning(), time.getMin(), time.getMax()
    time.setBinning(knotbinning, '%s_knotbinning' % pfx)
    time.setBinning(oldbinning)
    time.setRange(lo, hi)
    del knotbinning
    del oldbinning
    del lo
    del hi
    # create the knot coefficients
    coefflist = RooArgList()
    i = 0
    for v in mycoeffs:
        if floatParams:
            coefflist.add(WS(ws, RooRealVar('%s_SplineAccCoeff%u' % (pfx, i),
                'v_{%u}' % (i+1), v, 0., 3.)))
        else:
            coefflist.add(WS(ws, RooConstVar('%s_SplineAccCoeff%u' % (pfx, i),
                'v_{%u}' % (i+1), v)))
        i = i + 1
    del mycoeffs
    coefflist.add(one)
    i = i + 1
    myknots.append(time.getMax())
    myknots.reverse()
    fudge = (myknots[0] - myknots[1]) / (myknots[2] - myknots[1])
    lastmycoeffs = RooArgList(
            WS(ws, RooConstVar('%s_SplineAccCoeff%u_coeff0' % (pfx, i),
                '%s_SplineAccCoeff%u_coeff0' % (pfx, i), 1. - fudge)),
            WS(ws, RooConstVar('%s_SplineAccCoeff%u_coeff1' % (pfx, i),
                '%s_SplineAccCoeff%u_coeff1' % (pfx, i), fudge)))
    del myknots
    coefflist.add(WS(ws, RooPolyVar(
        '%s_SplineAccCoeff%u' % (pfx, i), 'v_{%u}' % (i+1),
        coefflist.at(coefflist.getSize() - 2), lastmycoeffs)))
    del i
    if debug:
        print 'DEBUG: Spline Coeffs: %s' % str([
            coefflist.at(i).getVal() for i in xrange(0, coefflist.getSize())
            ])
    # create the spline itself
    tacc = WS(ws, RooCubicSplineFun('%s_SplineAcceptance' % pfx,
        '%s_SplineAcceptance' % pfx, time, '%s_knotbinning' % pfx,
        coefflist))
    del lastmycoeffs
    if not floatParams:
        # make sure the acceptance is <= 1 for generation
        m = max([coefflist.at(j).getVal() for j in
            xrange(0, coefflist.getSize())])
#.........这里部分代码省略.........