diff --git a/python/sandbox.py b/python/sandbox.py
index 7d15581..f0db172 100644
--- a/python/sandbox.py
+++ b/python/sandbox.py
@@ -8,7 +8,10 @@ import config
 import microphone
 import dsp
 import led
+import melbank
 from scipy.ndimage.filters import gaussian_filter1d
+from scipy.signal import argrelextrema
+
 
 
 def rainbow(length, speed=1.0 / 5.0):
@@ -58,9 +61,9 @@ def rainbow_gen(length, speed=1./5., center=0.5, width=0.5, f=[1, 1, 1]):
     dt = 2.0 * np.pi / length
     t = time.time() * speed
     phi = 2.0 / 3.0 * np.pi
-    def r(t): return np.clip(np.sin(f[0] * t + 0. * phi) * width + center, 0., 1.)
-    def g(t): return np.clip(np.sin(f[1] * t + 1. * phi) * width + center, 0., 1.)
-    def b(t): return np.clip(np.sin(f[2] * t + 2. * phi) * width + center, 0., 1.)
+    def r(t): return np.clip(np.sin(f[0] * t + 1. * phi) * width + center, 0., 1.)
+    def g(t): return np.clip(np.sin(f[1] * t + 2. * phi) * width + center, 0., 1.)
+    def b(t): return np.clip(np.sin(f[2] * t + 3. * phi) * width + center, 0., 1.)
     x = np.tile(0.0, (length, 3))
     for i in range(length):
         x[i][0] = r(i * dt + t)
@@ -72,7 +75,7 @@ def rainbow_gen(length, speed=1./5., center=0.5, width=0.5, f=[1, 1, 1]):
 _time_prev = time.time() * 1000.0
 """The previous time that the frames_per_second() function was called"""
 
-_fps = dsp.ExponentialFilter(val=config.FPS, alpha_decay=0.05, alpha_rise=0.05)
+_fps = dsp.ExpFilter(val=config.FPS, alpha_decay=0.05, alpha_rise=0.05)
 """The low-pass filter used to estimate frames-per-second"""
 
 
@@ -123,7 +126,7 @@ def update_plot_1(x, y):
     global curves, p1
     colors = rainbow(config.N_CURVES) * 255.0
     for i in range(config.N_CURVES):
-        curves[i].setPen((colors[i][0], colors[i][1], colors[i][2]))
+        #curves[i].setPen((colors[i][0], colors[i][1], colors[i][2]))
         curves[i].setData(x=x, y=y[i])
     p1.autoRange()
     p1.setRange(yRange=(0.0, 2.0))
@@ -170,14 +173,14 @@ def leak_saturated_pixels(pixels):
     return pixels
 
 
-_EA_norm = dsp.ExponentialFilter(np.tile(1e-4, config.N_PIXELS), 0.01, 0.25)
+_EA_norm = dsp.ExpFilter(np.tile(1e-4, config.N_PIXELS), 0.01, 0.85)
 """Onset energy per-bin normalization constants
 
 This filter is responsible for individually normalizing the onset bin energies.
 This is used to obtain per-bin automatic gain control.
 """
 
-_EA_smooth = dsp.ExponentialFilter(np.tile(1.0, config.N_PIXELS), 0.25, 0.80)
+_EA_smooth = dsp.ExpFilter(np.tile(1.0, config.N_PIXELS), 0.25, 0.80)
 """Asymmetric exponential low-pass filtered onset energies
 
 This filter is responsible for smoothing the displayed onset energies.
@@ -216,8 +219,8 @@ def update_leds_6(y):
     return pixels
 
 
-_EF_norm = dsp.ExponentialFilter(np.tile(1.0, config.N_PIXELS), 0.05, 0.9)
-_EF_smooth = dsp.ExponentialFilter(np.tile(1.0, config.N_PIXELS), 0.08, 0.9)
+_EF_norm = dsp.ExpFilter(np.tile(1.0, config.N_PIXELS), 0.05, 0.9)
+_EF_smooth = dsp.ExpFilter(np.tile(1.0, config.N_PIXELS), 0.08, 0.9)
 _prev_energy = 0.0
 
 
@@ -236,8 +239,25 @@ def update_leds_5(y):
     return pixels
 
 
-_energy_norm = dsp.ExponentialFilter(10.0, alpha_decay=.15, alpha_rise=.9)
-_energy_smooth = dsp.ExponentialFilter(10.0, alpha_decay=0.1, alpha_rise=0.8)
+_prev_E = np.tile(.1, config.N_PIXELS)
+_EF_N = dsp.ExpFilter(np.tile(.1, config.N_PIXELS), 0.01, 0.9)
+def update_leds_5(y):
+    global _prev_E
+    y = np.copy(y)
+    current_E = y**2.0
+    EF = current_E - _prev_E
+    _prev_E = np.copy(current_E)
+    EF[EF < 0.0] = 0.0
+    _EF_N.update(EF)
+    EF /= _EF_N.value
+    EF[current_E < 0.02] = 0.0
+    return EF
+
+
+
+
+_energy_norm = dsp.ExpFilter(10.0, alpha_decay=.15, alpha_rise=.9)
+_energy_smooth = dsp.ExpFilter(10.0, alpha_decay=0.1, alpha_rise=0.8)
 
 
 # Modulate brightness by relative average rectified onset flux
@@ -284,9 +304,11 @@ def update_leds_2(y):
 
 def update_leds_1(y):
     """Display the raw onset spectrum on the LED strip"""
-    return np.copy(y)**0.5
+    return np.copy(y)**0.75
 
 
+YS_peak = dsp.ExpFilter(1.0, alpha_decay=0.005, alpha_rise=0.95)
+
 def microphone_update(stream):
     global y_roll
     # Retrieve new audio samples and construct the rolling window
@@ -295,39 +317,52 @@ def microphone_update(stream):
     y_roll = np.roll(y_roll, -1, axis=0)
     y_roll[-1, :] = np.copy(y)
     y_data = np.concatenate(y_roll, axis=0)
-    # Calculate onset detection functions
-    SF, NWPD, RCD = dsp.onset(y_data)
-    # Apply Gaussian blur to improve agreement between onset functions
-    SF = gaussian_filter1d(SF, 1.0)
-    NWPD = gaussian_filter1d(NWPD, 1.0)
-    RCD = gaussian_filter1d(RCD, 1.0)
-    # Update and normalize peak followers
-    SF_peak.update(np.max(SF))
-    NWPD_peak.update(np.max(NWPD))
-    RCD_peak.update(np.max(RCD))
-    SF /= SF_peak.value
-    NWPD /= NWPD_peak.value
-    RCD /= RCD_peak.value
-    # Normalize and update onset spectrum
-    # onset = np.sqrt(SF**2.0 + NWPD**2.0 + RCD**2.0)
-    # onset = SF * NWPD * RCD
-    onset = SF + NWPD + RCD
-    # onset = SF + RCD
-    onset_peak.update(np.max(onset))
-    onset /= onset_peak.value
-    onsets.update(onset)
-    # Map the onset values to LED strip pixels
-    if len(onsets.value) != config.N_PIXELS:
-        onset_values = interpolate(onsets.value, config.N_PIXELS)
-    else:
-        onset_values = np.copy(onsets.value)
-    brightness = led_visualization(onset_values)
+    # # Calculate onset detection functions
+    # SF, NWPD, RCD = dsp.onset(y_data)
+    # # Apply Gaussian blur to improve agreement between onset functions
+    # SF = gaussian_filter1d(SF, 1.0)
+    # NWPD = gaussian_filter1d(NWPD, 1.0)
+    # RCD = gaussian_filter1d(RCD, 1.0)
+    # # Update and normalize peak followers
+    # SF_peak.update(np.max(SF))
+    # NWPD_peak.update(np.max(NWPD))
+    # RCD_peak.update(np.max(RCD))
+    # SF /= SF_peak.value
+    # NWPD /= NWPD_peak.value
+    # RCD /= RCD_peak.value
+    # # Normalize and update onset spectrum
+    # onset = SF + NWPD + RCD
+    # onset_peak.update(np.max(onset))
+    # onset /= onset_peak.value
+    # onsets.update(onset)
+    # # Map the onset values to LED strip pixels
+    # if len(onsets.value) != config.N_PIXELS:
+    #     onset_values = interpolate(onsets.value, config.N_PIXELS)
+    # else:
+    #     onset_values = np.copy(onsets.value)
+    # brightness = led_visualization(onset_values)
+
+    XS, YS = dsp.fft(y_data, window=np.hamming)
+    YS = YS[XS >= 0.0]
+    XS = XS[XS >= 0.0]
+    YS = np.atleast_2d(np.abs(YS)).T * dsp.mel_y.T
+    YS = np.sum(YS, axis=0)**2.0
+    #YS = gaussian_filter1d(YS, 2.0)
+    YS = np.diff(YS, n=0)
+    YS_peak.update(np.max(YS))
+    YS /= YS_peak.value
+
+    if len(YS) != config.N_PIXELS:
+        YS = interpolate(YS, config.N_PIXELS)
+
+    #YS = led_visualization(YS)
+    YS = led_vis3(YS)
+
     # Plot the onsets
-    plot_x = np.array(range(1, len(onsets.value) + 1))
-    plot_y = [0*onsets.value**i for i in np.linspace(2.0, 0.25, config.N_CURVES)]
-    if brightness is not None:
-        plot_y = np.array([brightness, onsets.value])
-    #plot_y = brightness
+    #plot_x = np.array(range(1, len(onsets.value) + 1))
+    plot_x = np.array(range(1, len(YS) + 1))
+    #plot_y = [onsets.value**i for i in np.linspace(2.0, 0.25, config.N_CURVES)]
+    plot_y = [YS]
     update_plot_1(plot_x, plot_y)
     app.processEvents()
     print('FPS {:.0f} / {:.0f}'.format(frames_per_second(), config.FPS))
@@ -336,7 +371,7 @@ def microphone_update(stream):
 # Create plot and window
 app = QtGui.QApplication([])
 win = pg.GraphicsWindow('Audio Visualization')
-win.resize(800, 600)
+win.resize(300, 200)
 win.setWindowTitle('Audio Visualization')
 # Create plot 1 containing config.N_CURVES
 p1 = win.addPlot(title='Onset Detection Function')
@@ -354,20 +389,19 @@ pixels = np.tile(0.0, config.N_PIXELS)
 color = rainbow(config.N_PIXELS) * 255.0
 
 # Tracks average onset spectral energy
-onset_energy = dsp.ExponentialFilter(1.0, alpha_decay=0.01, alpha_rise=0.65)
+onset_energy = dsp.ExpFilter(1.0, alpha_decay=0.01, alpha_rise=0.65)
 
 # Tracks the location of the spectral median
-median = dsp.ExponentialFilter(val=config.N_SUBBANDS / 2.0,
+median = dsp.ExpFilter(val=config.N_SUBBANDS / 2.0,
                                alpha_decay=0.1, alpha_rise=0.1)
 # Smooths the decay of the onset detection function
-onsets = dsp.ExponentialFilter(val=np.tile(0.0, (config.N_SUBBANDS)),
+onsets = dsp.ExpFilter(val=np.tile(0.0, (config.N_SUBBANDS)),
                                alpha_decay=0.15, alpha_rise=0.75)
-
 # Peak followers used for normalization
-SF_peak = dsp.ExponentialFilter(1.0, alpha_decay=0.01, alpha_rise=0.99)
-NWPD_peak = dsp.ExponentialFilter(1.0, alpha_decay=0.01, alpha_rise=0.99)
-RCD_peak = dsp.ExponentialFilter(1.0, alpha_decay=0.01, alpha_rise=0.99)
-onset_peak = dsp.ExponentialFilter(0.1, alpha_decay=0.002, alpha_rise=0.5)
+SF_peak = dsp.ExpFilter(1.0, alpha_decay=0.01, alpha_rise=0.99)
+NWPD_peak = dsp.ExpFilter(1.0, alpha_decay=0.01, alpha_rise=0.99)
+RCD_peak = dsp.ExpFilter(1.0, alpha_decay=0.01, alpha_rise=0.99)
+onset_peak = dsp.ExpFilter(0.1, alpha_decay=0.002, alpha_rise=0.5)
 
 # Number of audio samples to read every time frame
 samples_per_frame = int(config.MIC_RATE / config.FPS)
@@ -382,30 +416,47 @@ y_roll = np.random.rand(config.N_ROLLING_HISTORY, samples_per_frame) / 100.0
 # update_leds_5 = energy flux normalized per-bin spectrum (GAMMA = True)
 # update_leds_6 = energy average normalized per-bin spectrum (GAMMA = True)
 
-
 # Low pass filter for the LEDs being output to the strip
-pixels_filt = dsp.ExponentialFilter(np.tile(0., (config.N_PIXELS, 3)), .2, .8)
+pixels_filt = dsp.ExpFilter(np.tile(0., (config.N_PIXELS, 3)), .14, .9)
 
 
 def hyperbolic_tan(x):
     return 1.0 - 2.0 / (np.exp(2.0 * x) + 1.0)
 
+
+def bloom_peaks(x, width=3, blur_factor=1.0):
+    peaks = argrelextrema(x, np.greater)[0]
+    y = x * 0.0
+    if len(peaks) == 0:
+        return y
+    for peak in peaks:
+        min_idx = max(peak - width, 0)
+        max_idx = min(peak + width, len(x) - 1)
+        for i in range(min_idx, max_idx):
+            y[i] = x[i]
+    y = gaussian_filter1d(y, blur_factor)
+    return y
+
+
 # This is the function responsible for updating LED values
 # Edit this function to change the visualization
 def led_visualization(onset_values):
     # Visualizations that we want to use (normalized to ~[0, 1])
-    pixels_A = update_leds_6(onset_values)
-    pixels_B = update_leds_4(onset_values)
+    #pixels_A = update_leds_6(onset_values)
+    #pixels_B = update_leds_4(onset_values)
     # Combine the effects by taking the product
-    brightness = pixels_A * pixels_B
-    brightness = gaussian_filter1d(brightness, 1.0)**1.5
-    brightness = hyperbolic_tan(brightness)
+    #brightness = pixels_A #* pixels_B
+    brightness = update_leds_6(onset_values**2.0)
+    brightness = gaussian_filter1d(brightness, 4.0)
+    #brightness = hyperbolic_tan(brightness)
+    brightness = bloom_peaks(brightness)**2.
     # Combine pixels with color map
     color = rainbow_gen(onset_values.shape[0],
                         speed=1.,
                         center=0.5,
                         width=0.5,
-                        f=[1.0, 1.0, 1.]) * 255.0
+                        f=[1.1, .5, .2]) * 255.0
+    color = np.tile(255.0, (config.N_PIXELS, 3))
     # color = rainbow(onset_values.shape[0]) * 255.0
     pixels = (brightness * color.T).T
     pixels = leak_saturated_pixels(pixels)
@@ -418,6 +469,58 @@ def led_visualization(onset_values):
     return brightness
 
 
+mean_energy = dsp.ExpFilter(0.1, alpha_decay=0.05, alpha_rise=0.05)
+
+def led_vis2(x):
+    energy = np.mean(x**.5)
+    mean_energy.update(energy)
+    energy = energy / mean_energy.value - 1.0
+    edge = np.exp(-10 * np.linspace(0, 1, len(x)))
+    edge = edge + edge[::-1]
+    edge *= max(energy, 0)
+    edge /= 2.0
+
+    x = gaussian_filter1d(x, 3.0)
+    x = update_leds_6(x)
+    red = bloom_peaks(x**1.0, width=1, blur_factor=1.5)
+    green = bloom_peaks(x**1.0, width=2, blur_factor=0.5)
+    blue = bloom_peaks(x**1.0, width=1, blur_factor=0.5)
+    # Set LEDs
+    color = np.tile(0.0, (3, config.N_PIXELS))
+    color[0, :] = 1.0*edge + red*1.0
+    color[1, :] = 1.2*edge + green*1.0
+    color[2, :] = 1.5*edge + blue*1.0
+    color = color.T * 255.0
+    pixels_filt.update(color)
+    led.pixels = np.round(pixels_filt.value).astype(int)
+    led.update()
+    return (color[:, 0] + color[:, 1] + color[:, 2]) / (3. * 255.0)
+
+
+N = 60
+E = []
+for i in range(0, N):
+    alpha_decay = 0.01 * (float(i + 1) / (N + 1.0))**2.0
+    alpha_rise = alpha_decay
+    E.append(dsp.ExpFilter(.1, alpha_decay, alpha_rise))
+
+def led_vis3(x):
+    energy = np.mean(x**.5)
+    pixels = np.tile(0.0, config.N_PIXELS)
+    for i in range(N):
+        E[i].update(energy)
+        pixels[i] = hyperbolic_tan(max(energy / E[i].value - 1.0, 0))
+
+    color = np.tile(0.0, (3, config.N_PIXELS))
+    color[0, :] = pixels
+    color[1, :] = pixels
+    color[2, :] = pixels
+    color = color.T * 255.0
+    pixels_filt.update(color)
+    led.pixels = np.round(pixels_filt.value).astype(int)
+    led.update()
+    return (color[:, 0] + color[:, 1] + color[:, 2]) / (3. * 255.0)
+
 if __name__ == '__main__':
     led.update()
     microphone.start_stream(microphone_update)