from ROOT import gRandom, TF2
from math import cos, sqrt

x_a = 0.01
x_b = 1
y_a = 0
y_b = 1
n = 1000
m = 10
N = 10000

def norm():
    return  1.0 / (2 * x_a ** 2) - 1.0 / (2 * x_b ** 2)

def func(x,y):
    return 1.0 / (x ** 3 * (y + cos(x) ** 2 / 3))

def g(x,y):
    return 1.0 / (x ** 3) / norm()

def F(x,y):
    return func(x,y) / g(x,y)

def integrate(func, xlim, ylim):
    xa, xb = xlim
    ya, yb = ylim
    x = gRandom.Uniform(xa, xb)
    y = gRandom.Uniform(ya, yb)
    return (xb - xa) * (yb - ya) * func(x,y)

def fndirect(pp):
    xa = pp[0]
    xb = pp[1]
    ya = pp[2]
    yb = pp[3]
    ua = 1.0 / (2 * xb ** 2)
    ub = 1.0 / (2 * xa ** 2)
    u = gRandom.Uniform(ub, ua)
    r1 = 1.0 / sqrt(2 * u + 1.0 / xa ** 2)
    r2 = gRandom.Uniform(ya, yb)
    return r1, r2

def meanmethod():
    res = []
    print "MEAN METHOD"
    for i in range(m):
        cres = []
        for j in range(n):
        r = (integrate(func, (x_a, x_b), (y_a, y_b)))
        res.append(r)
        cres.append(r)
        avg = sum(res) / len(res)
        var = (sum(r ** 2 for r in res) - avg ** 2) / len(res)
        cavg = sum(cres) / len(cres)
        cvar = (sum(r ** 2 for r in cres) - cavg ** 2) / len(cres)
        rms = sqrt(var / len(res))
        crms = sqrt(cvar / len(cres))
        print "%d iteration: average is %e, rms is %e " % ( i, avg, rms)
        print "%d iteration: caverage is %e, crms is %e " % ( i, cavg, crms)

def sampling():
    res = []
    print "SAMPLING METHOD"
    for i in range(m):
    cres = []
    for j in range(n):
        r1, r2 = fndirect([x_a, x_b, y_a, y_b])
        r = F(r1, r2)
            res.append(r)
        cres.append(r)
    avg = sum(res) / len(res)
        var = (sum(r ** 2 for r in res) - avg ** 2) / len(res)
        cavg = sum(cres) / len(cres)
        cvar = (sum(r ** 2 for r in cres) - cavg ** 2) / len(cres)
        rms = sqrt(var / len(res))
        crms = sqrt(cvar / len(cres))
        print "%d iteration: average is %e, rms is %e " % ( i, avg, rms)
        print "%d iteration: caverage is %e, crms is %e " % ( i, cavg, crms)

def drawfunc( xx, pp):
    return func(xx[0], xx[1])
def drawF(xx,pp):
    return F(xx[0], xx[1])
meanmethod()
sampling()
h = TF2("f", drawfunc, x_a, x_b, y_a, y_b, 1)
h.Draw("surf")
raw_input()
h1 = TF2("F", drawF, x_a, x_b, y_a, y_b, 1)
h1.Draw("surf")
raw_input()