# Path: fake_spectra/fluxstatistics.py
# def obs_mean_tau(redshift):
# def mean_flux(tau, mean_flux_desired, tol = 1e-5, thresh=1e30):
# def flux_pdf(tau, nbins=20, mean_flux_desired=None):
# def _powerspectrum(inarray, axis=-1):
# def _window_function(k, *, R, dv):
# def flux_power(tau, vmax, spec_res = 8, mean_flux_desired=None, window=True):
# def _flux_power_bins(vmax, npix):


import math
import numpy as np
from fake_spectra import fluxstatistics as stat
    assert np.all(np.abs(hist - np.array([ 16. ,   2.2,   0.6,   0.2,   0.2,   0.2,   0.2,   0. ,   0. , 0. ,   0.2,   0. ,   0. ,   0. ,   0. ,   0. ,   0. ,   0. , 0. ,   0.2])) < 1e-3)

def testPowerspectrum():
    """Test the core FFT routine"""
    xx = np.tile(np.arange(0,1,0.01)**2,(10,1))
#     xx = np.arange(0,1,0.01)**2
    assert np.shape(xx)[0] == 10
    fpk = stat._powerspectrum(np.exp(-xx), axis=1) * (2 * math.pi)
    #Check that Parseval's theorem is true, accounting for the negative frequency modes not included in the DFT.
    for ff in fpk:
        dpower = (np.sum(ff)+np.sum(ff[1:]))
        assert abs(dpower - 2 * math.pi *np.sum(np.exp(-xx[0,:])**2)/np.shape(xx)[1]) < 0.1
    xx = np.linspace(1,51,200)
    inn = np.sin(2*math.pi*xx)
    #This will be a delta function centered at 50
    ff = stat._powerspectrum(inn)
    assert np.where(np.max(ff) == ff)[0][0] == 50

def testFluxPower():
    """Check that taking the mean of the flux power spectra is working correctly"""
    #Check both even and odd binning
    for bb in (200, 201):
        xx = np.linspace(0,51,bb)
        inn = np.sin(2*math.pi*xx)+1.5
        ff = stat._powerspectrum(np.exp(-inn)-1)
        #Construct some optical depths offset from each other
        taus = np.zeros((9,bb))
        taus = np.vstack([inn,]*10)
        bins, power = stat.flux_power(taus, vmax=1.,spec_res = 0.01)
        power/=12.5569

# Path: fake_spectra/fluxstatistics.py
# def obs_mean_tau(redshift):
# def mean_flux(tau, mean_flux_desired, tol = 1e-5, thresh=1e30):
# def flux_pdf(tau, nbins=20, mean_flux_desired=None):
# def _powerspectrum(inarray, axis=-1):
# def _window_function(k, *, R, dv):
# def flux_power(tau, vmax, spec_res = 8, mean_flux_desired=None, window=True):
# def _flux_power_bins(vmax, npix):


# Path: fake_spectra/tests/test_statistics.py
import math
import numpy as np
from fake_spectra import fluxstatistics as stat

    assert np.all(np.abs(hist - np.array([ 16. ,   2.2,   0.6,   0.2,   0.2,   0.2,   0.2,   0. ,   0. , 0. ,   0.2,   0. ,   0. ,   0. ,   0. ,   0. ,   0. ,   0. , 0. ,   0.2])) < 1e-3)

def testPowerspectrum():
    """Test the core FFT routine"""
    xx = np.tile(np.arange(0,1,0.01)**2,(10,1))
#     xx = np.arange(0,1,0.01)**2
    assert np.shape(xx)[0] == 10
    fpk = stat._powerspectrum(np.exp(-xx), axis=1) * (2 * math.pi)
    #Check that Parseval's theorem is true, accounting for the negative frequency modes not included in the DFT.
    for ff in fpk:
        dpower = (np.sum(ff)+np.sum(ff[1:]))
        assert abs(dpower - 2 * math.pi *np.sum(np.exp(-xx[0,:])**2)/np.shape(xx)[1]) < 0.1
    xx = np.linspace(1,51,200)
    inn = np.sin(2*math.pi*xx)
    #This will be a delta function centered at 50
    ff = stat._powerspectrum(inn)
    assert np.where(np.max(ff) == ff)[0][0] == 50

def testFluxPower():
    """Check that taking the mean of the flux power spectra is working correctly"""
    #Check both even and odd binning
    for bb in (200, 201):
        xx = np.linspace(0,51,bb)
        inn = np.sin(2*math.pi*xx)+1.5
        ff = stat._powerspectrum(np.exp(-inn)-1)
        #Construct some optical depths offset from each other
        taus = np.zeros((9,bb))
        taus = np.vstack([inn,]*10)
        bins, power = stat.flux_power(taus, vmax=1.,spec_res = 0.01)
        power/=12.5569

        wind = stat._window_function(bins[1:],R=0.01,dv=1/np.size(xx))