import sys
import numpy as np
from numpy import sin, cos, tan, pi
from pprint import pprint


def SmoothingBySimpleAverage(y, n):
    n2 = int(n / 2);
    ndata = len(y);
    ys = []
    for i in range(0, ndata):
        c = 0;
        ys.append(0.0);
        for k in range(i - n2, i + n2 + 1):
            if k < 0 or k >= ndata:
                continue
            ys[i] += y[k]
            c += 1
        if c > 0:
            ys[i] /= c;
        else:
            ys[i] = y[i]
    return ys;

def SmoothingByPolynomialFit(y, n):
    m = int(n / 2);
    n = m + m + 1
    W23 = (4.0 * m * m - 1.0) * (2.0 * m + 3.0) / 3.0;
    w23j = [0.0]*n
    for j in range(-m, m+1):
        w23j[j + m] = (3.0 * m * (m+1.0) - 1.0 - 5.0 * j * j) / W23

    ndata = len(y)
    ys = []
    for i in range(0, ndata):
        c = 0.0;
        ys.append(0.0);
        for j in range(-m, m+1):
            k = i + j
            if k < 0 or k >= ndata:
                continue
            ys[i] += w23j[j+m] * y[k]
            c += w23j[j+m]
        if c > 0:
            ys[i] /= c
        else:
            ys[i] = y[i]
    return ys;

def polynomial_fit_by_range(y, ny, nsmooth, i0):
    m = int(nsmooth / 2);
# ny = 10  nsmooth = 3  i0 = 2 => m = 2
#                       i0 = 7 => ny - 1 - i0 = 2 => m = 2
    if i0 < nsmooth:
        m = i0
    if ny - 1 - i0 < nsmooth:
        m = ny - 1 - i0
    nsmooth = m + m + 1

    W23 = (4.0 * m * m - 1.0) * (2.0 * m + 3.0) / 3.0;
    w23j = [0.0]*nsmooth
    for j in range(-m, m+1):
        w23j[j + m] = (3.0 * m * (m+1.0) - 1.0 - 5.0 * j * j) / W23

    ys = 0.0
    c  = 0.0
    for j in range(-m, m+1):
        ys  += w23j[j+m] * y[i0+j]
        c += w23j[j+m]

    ys /= c

    return ys

def polynomial_fit(y, nsmooth):
    ny = len(y)
    if ny <= 1:
        return y

    m = int(nsmooth / 2);
    nsmooth = m + m +1
    ys = []
    for i in range(ny):
        s = polynomial_fit_by_range(y, ny, nsmooth, i)
        ys.append(s)
    return ys

def simple_average(y, nsmooth):
    m = int(nsmooth / 2);
    nsmooth = m + m +1
    ys = []
    ny = len(y)
    for i in range(ny):
        if i == 0:
            ys.append(y[i])
        elif i == ny - 1:
            ys.append(y[ny - 1])
        elif i < nsmooth:
            ms = i
            ns = ms + ms + 1
            s = sum(y[0:ns]) / ns
            ys.append(s)
        elif i > ny - nsmooth - 1:
            ms = ny - 1 - i
            ns = ms + ms + 1
            s = sum(y[i - ms:i + ms + 1]) / ns
            ys.append(s)
        else:
            s = sum(y[i - m:i + m + 1]) / nsmooth
            ys.append(s)
    return ys

