# Path: lmatools/flash_stats.py
# def raw_moments(a, n_moments=5):
#     """ Given data in array a, return the raw moments,
#         sum(a^m) where m is 0...(n_moments-1).
#         Return is an array of shape (n_moments,).
#     """
#     sample_raw_moments = np.zeros(n_moments, dtype='f8')
#     for i in range(n_moments):
#         sample_raw_moments[i] = (a**i).sum()
#     sample_raw_moments[1:] /= sample_raw_moments[0]    
#     return sample_raw_moments
# 
# @coroutine        
# def raw_moments_for_parameter(parameter, preprocess=None, n_moments=5, output_target=None):
#     """ This coroutine builds up raw moments by streaming samples into the coroutine.
#         When the samples are exhausted by a GeneratorExit, An array of size n_moments 
#         will be sent to the output_target, and will contain the zeroth and n_moments-1 
#         higher moments. The higher moments are already divided by the zeroth moment.
#         
#         Receives a data array, and calculate basic statistics from the distribution
#         of the named parameter in that array. Optionally call the preprocess function 
#         on the parameter before calculating the moment.
#         
#     """
#     sample_raw_moments = np.zeros(n_moments, dtype='f8')
#     try:
#         while True:
#             data = (yield)
#             a = data[parameter]
#             if preprocess is not None:
#                 a = preprocess(a)
#             # calculate the sample moments
#             for i in range(n_moments):
#                 sample_raw_moments[i] += (a**i).sum()
#             
#     except GeneratorExit:
#         sample_raw_moments[1:] /= sample_raw_moments[0]
#         if output_target is not None:
#             output_target.send(sample_raw_moments)
# 
# def central_moments_from_raw(raw):
#     """ Returns ctr, std where
#             ctr = the zero through fourth central moments
#             std = (number, mean, variance, skewness, kurtosis)
#             
#             the zeroth and first central moments are set to zero.
#     """
#     ctr = np.zeros_like(raw)
#     # ctr[0] = 0
#     # ctr[1] = raw[1] - raw[1]
#     ctr[2] = raw[2] - raw[1]*raw[1]
#     ctr[3] = 2*(raw[1]**3) - 3*raw[1]*raw[2] + raw[3]
#     ctr[4] = -3*(raw[1]**4) + 6*raw[1]*raw[1]*raw[2] - 4*raw[3]*raw[1] + raw[4]
#     # number, mean, variance, skewness, kurtosis
#     std = raw[0], raw[1], ctr[2], ctr[3]/(ctr[2]**1.5), (ctr[4]/(ctr[2]*ctr[2]) - 3)
#     return ctr, std
# 
# @coroutine
# def events_flashes_receiver(target=None):
#     """ Passes along only flashes """
#     while True:
#         events, flashes = (yield)
#         if target is not None:
#             target.send(flashes)
# 
# Path: lmatools/grid/make_grids.py
# def time_edges(start_time, end_time, frame_interval):
#     """ Return lists cooresponding the start and end times of frames lasting frame_interval
#         between start_time and end_time. The last interval may extend beyond end_time, but
#         by no more than frame_interval. This makes each frame the same length.
# 
#         returns t_edges, duration, where t_edges is a list of datetime objects, and
#         duration is the total duration between the start and end times (and not the duration
#         of all frames)
#     """
#     frame_dt = timedelta(0, frame_interval, 0)
#     duration = end_time - start_time
#     n_frames = int(np.ceil(to_seconds(duration) / to_seconds(frame_dt)))
#     t_edges = [start_time + i*frame_dt for i in range(n_frames+1)]
#     return t_edges, duration
# 
# def seconds_since_start_of_day(start_time, t):
#     """ For each datetime object t, return the number of seconds elapsed since the
#         start of the date given by start_time. Only the date part of start_time is used.
#     """
#     ref_date = start_time.date()
#     t_ref = datetime(ref_date.year, ref_date.month, ref_date.day)
#     t_edges_seconds = [to_seconds(edge - t_ref) for edge in t]
#     return t_ref, t_edges_seconds


import os
import numpy as np
import matplotlib.pyplot as plt
from stormdrain.pipeline import coroutine
from stormdrain.bounds import Bounds, BoundsFilter
from stormdrain.support.matplotlib.formatters import SecDayFormatter
from matplotlib.ticker import MultipleLocator
from lmatools.flash_stats import raw_moments, raw_moments_for_parameter, central_moments_from_raw, events_flashes_receiver
from lmatools.grid.make_grids import time_edges, seconds_since_start_of_day
from stormdrain.support.matplotlib.poly_lasso import LassoFilter
            lines.sort()
            for l in lines:
                print(l)
            fig = self.plot()
            fig.savefig(self.outdir+'/moment-energy-timeseries.pdf')
        

def flash_size_stats_for_intervals(t_start, t_end, dt, lat_bounds=None, lon_bounds=None, outdir=None):
    """ Create a pipeline to process an arbitrary number of flashes which are 
        accumulated in dt-length intervals between t_start and t_end.
        
        t_start, t_end: datetime.datetime instances
        dt: datetime.timedelta
        
        Returns:
        printer_of_stats: 
        all_frame_targets: a list of targets that receive (events, flashes) from the flash file reader

        all_frame_targets can be passed as other_analysis_targets to lmatools.flash_stats.plot_spectra_for_files.
        
        When all data have been read (i.e., plot_spectra_for_files has returned), the accumulated statistics are
        retrieved by closing the targets returned by this function, i.e.
            for target in timeframe_targets:
                target.close()
            printer_of_stats.close()
        
        plot_spectra_for_files itself could be refactored to make use of these same time windows natively, but
        that's a project for later. The dependency here on stormdrain is another issue; it may make sense to have
        lmatools depend on stormdrain.
    """

# Path: lmatools/flash_stats.py
# def raw_moments(a, n_moments=5):
#     """ Given data in array a, return the raw moments,
#         sum(a^m) where m is 0...(n_moments-1).
#         Return is an array of shape (n_moments,).
#     """
#     sample_raw_moments = np.zeros(n_moments, dtype='f8')
#     for i in range(n_moments):
#         sample_raw_moments[i] = (a**i).sum()
#     sample_raw_moments[1:] /= sample_raw_moments[0]    
#     return sample_raw_moments
# 
# @coroutine        
# def raw_moments_for_parameter(parameter, preprocess=None, n_moments=5, output_target=None):
#     """ This coroutine builds up raw moments by streaming samples into the coroutine.
#         When the samples are exhausted by a GeneratorExit, An array of size n_moments 
#         will be sent to the output_target, and will contain the zeroth and n_moments-1 
#         higher moments. The higher moments are already divided by the zeroth moment.
#         
#         Receives a data array, and calculate basic statistics from the distribution
#         of the named parameter in that array. Optionally call the preprocess function 
#         on the parameter before calculating the moment.
#         
#     """
#     sample_raw_moments = np.zeros(n_moments, dtype='f8')
#     try:
#         while True:
#             data = (yield)
#             a = data[parameter]
#             if preprocess is not None:
#                 a = preprocess(a)
#             # calculate the sample moments
#             for i in range(n_moments):
#                 sample_raw_moments[i] += (a**i).sum()
#             
#     except GeneratorExit:
#         sample_raw_moments[1:] /= sample_raw_moments[0]
#         if output_target is not None:
#             output_target.send(sample_raw_moments)
# 
# def central_moments_from_raw(raw):
#     """ Returns ctr, std where
#             ctr = the zero through fourth central moments
#             std = (number, mean, variance, skewness, kurtosis)
#             
#             the zeroth and first central moments are set to zero.
#     """
#     ctr = np.zeros_like(raw)
#     # ctr[0] = 0
#     # ctr[1] = raw[1] - raw[1]
#     ctr[2] = raw[2] - raw[1]*raw[1]
#     ctr[3] = 2*(raw[1]**3) - 3*raw[1]*raw[2] + raw[3]
#     ctr[4] = -3*(raw[1]**4) + 6*raw[1]*raw[1]*raw[2] - 4*raw[3]*raw[1] + raw[4]
#     # number, mean, variance, skewness, kurtosis
#     std = raw[0], raw[1], ctr[2], ctr[3]/(ctr[2]**1.5), (ctr[4]/(ctr[2]*ctr[2]) - 3)
#     return ctr, std
# 
# @coroutine
# def events_flashes_receiver(target=None):
#     """ Passes along only flashes """
#     while True:
#         events, flashes = (yield)
#         if target is not None:
#             target.send(flashes)
# 
# Path: lmatools/grid/make_grids.py
# def time_edges(start_time, end_time, frame_interval):
#     """ Return lists cooresponding the start and end times of frames lasting frame_interval
#         between start_time and end_time. The last interval may extend beyond end_time, but
#         by no more than frame_interval. This makes each frame the same length.
# 
#         returns t_edges, duration, where t_edges is a list of datetime objects, and
#         duration is the total duration between the start and end times (and not the duration
#         of all frames)
#     """
#     frame_dt = timedelta(0, frame_interval, 0)
#     duration = end_time - start_time
#     n_frames = int(np.ceil(to_seconds(duration) / to_seconds(frame_dt)))
#     t_edges = [start_time + i*frame_dt for i in range(n_frames+1)]
#     return t_edges, duration
# 
# def seconds_since_start_of_day(start_time, t):
#     """ For each datetime object t, return the number of seconds elapsed since the
#         start of the date given by start_time. Only the date part of start_time is used.
#     """
#     ref_date = start_time.date()
#     t_ref = datetime(ref_date.year, ref_date.month, ref_date.day)
#     t_edges_seconds = [to_seconds(edge - t_ref) for edge in t]
#     return t_ref, t_edges_seconds


# Path: lmatools/lasso/energy_stats.py
import os
import numpy as np
import matplotlib.pyplot as plt
from stormdrain.pipeline import coroutine
from stormdrain.bounds import Bounds, BoundsFilter
from stormdrain.support.matplotlib.formatters import SecDayFormatter
from matplotlib.ticker import MultipleLocator
from lmatools.flash_stats import raw_moments, raw_moments_for_parameter, central_moments_from_raw, events_flashes_receiver
from lmatools.grid.make_grids import time_edges, seconds_since_start_of_day
from stormdrain.support.matplotlib.poly_lasso import LassoFilter

            lines.sort()
            for l in lines:
                print(l)
            fig = self.plot()
            fig.savefig(self.outdir+'/moment-energy-timeseries.pdf')
        

def flash_size_stats_for_intervals(t_start, t_end, dt, lat_bounds=None, lon_bounds=None, outdir=None):
    """ Create a pipeline to process an arbitrary number of flashes which are 
        accumulated in dt-length intervals between t_start and t_end.
        
        t_start, t_end: datetime.datetime instances
        dt: datetime.timedelta
        
        Returns:
        printer_of_stats: 
        all_frame_targets: a list of targets that receive (events, flashes) from the flash file reader

        all_frame_targets can be passed as other_analysis_targets to lmatools.flash_stats.plot_spectra_for_files.
        
        When all data have been read (i.e., plot_spectra_for_files has returned), the accumulated statistics are
        retrieved by closing the targets returned by this function, i.e.
            for target in timeframe_targets:
                target.close()
            printer_of_stats.close()
        
        plot_spectra_for_files itself could be refactored to make use of these same time windows natively, but
        that's a project for later. The dependency here on stormdrain is another issue; it may make sense to have
        lmatools depend on stormdrain.
    """

    t_edges, duration = time_edges(t_start, t_end, dt.total_seconds())