obs

#!/usr/bin/env python

import numpy as np
import time
import math
import sys
from pymlab import config

cfg = config.Config(
    i2c = {
        'port': 0
    },
    bus = [
        {
            "type": "i2chub",
            "address": 0x70,
            "children": [
                { "name":"clkgen", "type":"clkgen01", "channel": 5, },
            ],
        },
    ],
)


def reduce(a, div=4, start=0, stop=2):
    return np.sum(np.sort(a, axis=0)[len(a)*start/div:len(a)*stop/div], axis=0) / (len(a) / 2)


def coroutine():
    def sleep(time):
        for i in xrange(int(math.ceil(float(time) / row_duration))):
            yield

    cfg.initialize()
    fgen = cfg.get_device('clkgen')
    fgen.reset()

    for i in sleep(3.0):
        yield

    fgen = cfg.get_device('clkgen')
    fgen.recall_nvm()

    for i in sleep(2.0):
        yield

    #freqs = [29200000, 29300000, 29400000, 29500000, 29600000, 29700000, 29800000]
    freqs = xrange(28000000, 31000000, 200000)

    nmeas_rows = int(math.ceil(float(1.0) / row_duration))
    arr = np.zeros(nmeas_rows, dtype=np.float32)

    avg_win = []

    while True:
#        avg = 0.0
        measurements = []

        for freq in freqs:
            fgen.recall_nvm()
            for i in sleep(0.2):
                yield
            freq_mhz = float(freq) / 1000000
            fgen.set_freq(10., freq_mhz * 2)
            for i in sleep(0.5):
                yield

            row, _s, _n = yield

            emit_event("mlab.aabb_event.measurement_area", (row, row + nmeas_rows, 0, 4096, "%f MHz" % (freq_mhz,)))

            for i in xrange(nmeas_rows):
                _r, _s, noise_lvl = yield
                arr[i] = noise_lvl

            m = len(arr) / 8
            #noise_lvl_sum = reduce(np.concatenate((arr[m:m*3], arr[m*5:m*7])))
            noise_lvl_sum = np.min(np.concatenate((arr[m:m*3], arr[m*5:m*7])))

            #sys.stdout.write("\t%f\t%f\t%f\n" % (time.time(), freq_mhz, noise_lvl_sum))
            #sys.stdout.flush()

            #avg += noise_lvl_sum
            measurements.append(noise_lvl_sum)


        #avg /= len(freqs)
        measurements.sort()
	val = sum(measurements[2:5])/3

        avg_win.append(val)
        avg_win = avg_win[-3:]

        sys.stdout.write("\t%f\t%f\n" % (time.time(), sum(avg_win)/3))
        sys.stdout.flush()


coroutine_inst = None

def run(row, spectrum):
    global coroutine_inst

    if coroutine_inst is None:
        coroutine_inst = coroutine()
        coroutine_inst.send(None)

    spectrum = np.log10(spectrum) * 10
    noise_lvl = reduce(spectrum)
    plot("noise", noise_lvl / 300)
    try:
        coroutine_inst.send((row, spectrum, noise_lvl))
    except StopIteration:
        pass

def process(sig_input, nbins, overlap):
    window = 0.5 * (1.0 - np.cos((2 * math.pi * np.arange(nbins)) / nbins))

    process_row = 0
    ringbuf = np.zeros(nbins * 4, dtype=np.complex64)
    ringbuf_edge = nbins
    readsize = nbins - overlap

    while True:
        if (ringbuf_edge + readsize > len(ringbuf)):
            ringbuf[0:overlap] = ringbuf[ringbuf_edge - overlap:ringbuf_edge]
            ringbuf_edge = overlap

        ringbuf[ringbuf_edge:ringbuf_edge + readsize] = sig_input.read(readsize)
        ringbuf_edge += readsize

        signal = ringbuf[ringbuf_edge - nbins:ringbuf_edge]

        spectrum = np.absolute(np.fft.fft(np.multiply(signal, window)))
        spectrum = np.concatenate((spectrum[nbins/2:nbins], spectrum[0:nbins/2]))

        run(process_row, spectrum)
        process_row = process_row + 1


class RawSigInput:
    def __init__(self, sample_rate, no_channels, dtype, file):
        self.sample_rate = sample_rate
        self.no_channels = no_channels
        self.dtype = dtype
        self.file = file

    def read(self, frames):
        read_len = frames * self.dtype.itemsize * self.no_channels
        string = ""

        while len(string) < read_len:
            string += self.file.read(read_len - len(string))

        if self.no_channels == 1:
            return np.fromstring(string, dtype=self.dtype).astype(np.float32)
        elif self.no_channels == 2 and self.dtype == np.dtype(np.float32):
            return np.fromstring(string, dtype=np.complex64)
        else:
            raise NotImplementedError("unimplemented no of channels and type combination")

    def start(self):
        pass

    def __str__(self):
        return "raw input from '%s'" % self.file.name

if __name__ == "__main__":
    global row_duration, plot, emit_event
    plot = lambda a, b: None
    emit_event = lambda a, b: None
    nbins = 4096
    overlap = 3072
    sig_input = RawSigInput(48000, 2, np.dtype(np.float32), sys.stdin)
    row_duration = float(nbins - overlap) / sig_input.sample_rate
    process(sig_input, 4096, 3072)