// TODO: declare audio ended so we can ignore any incomplete bass frequencies

sclass ContinuousOscillators_v2 extends Meta {

  macro dprint { if (debug) printVars } bool debug;

  settable bool bassFirst;

  settable bool printSubtractions;

  int channels = 1;

  settable double sampleRate = 48000;

  double currentSample;

  IAudioSample audio;              // integrated audio sample

  SubtractedAudio subtractedAudio; // audio with deleted frequencies

  settable bool keepRendered;

  transient double[] rendered;

  bool renderUsingSine; // or square wave

  bool realOnly, imagOnly; // for testing

  double maxRelFreq = 0.5; // maximum oscillator frequency relative to sample rate

  event escalating(Oscillator o);

  long haarFeaturesChecked;

  persistable class Oscillator extends HasKey<Frequency> {

    double currentPeriodStart;

    // all recorded intensities

    new TreeMap<DoubleRange, Channels<Complex>> intensities;

    double interval; // length of period in samples

    *(Frequency f) {

      super(f);

      interval = frequency().interval()*sampleRate;

    }

    Frequency frequency() { ret key(); }

    double interval() { ret interval; }

    scaffolded void record_impl(DoubleRange r, Channels<Complex> c) {

      intensities.put(r, c);

      subtractMeasured(r, c);

      currentPeriodStart = phase360();

      if (!bassFirst) escalate();

    }

    void escalate {

      escalating(this);

      rcallF(oscillators.higher(this), o -> { o?.process(currentPeriodStart); });

    }

    void subtractMeasured(DoubleRange r, Channels<Complex> c) {

      renderer from method;

      if (printSubtractions) {

        print("subtract " + r + " = " + c);

        enableScaffolding(renderer);

      }

      renderer.subtractPeriodFrom(subtractedAudio, r, frequency(), c.getSingleton());

      renderer().addPeriodTo(rendered, r, frequency(), c.getSingleton());

    }

    scaffolded void process(double upTo) {

      //printVars(+phase360(), +upTo);

      while (phase360() <= upTo)

        measure();

      subtractedAudioSanityCheck();

      escalate();

    }

    scaffolded void measure() {

      DoubleRange p1 = doubleRangeWithLength(currentPeriodStart, interval);

      AudioHaarFeature haar = new(relevantAudio(), p1);

      haarFeaturesChecked += 2;

      record(this, p1, haar.getComplex());

    }

    double phase360() { ret currentPeriodStart + interval(); }

  } // end of Oscillator

  IAudioSample relevantAudio() { ret subtractedAudio; }

  TreeHasKeyMap<Frequency, Oscillator> oscillators;

  *() {}

  *(IAudioSample audio, double *currentSample) { setAudio(audio); }

  *(IAudioSample audio) { setAudio(audio); } // start at time 0

  void setFrequencies(Iterable<Frequency> frequencies) {

    oscillators = new TreeHasKeyMap<Frequency, Oscillator>(map(f -> new Oscillator(f), filter(frequencies,

      f -> between(f!, 1, sampleRate*maxRelFreq))));

  }

  void setAudio(IAudioSample audio) {

    this.audio = audio;

    subtractedAudio = new SubtractedAudio(audio);

    if (keepRendered) rendered = new double[iceil(audio.length())];

    assertEquals("bounds", audio.bounds(), subtractedAudio.bounds());

  }

  int nFrequencies() { ret oscillators.size(); }

  scaffolded void stepTo(DoubleRange timeRange) {

    var o = first(oscillators);

    o?.process(timeRange.end);

  }

  void record(Oscillator o, DoubleRange r, Channels<Complex> c) {

    if (o == null) ret;

    o.record_impl(r, c);

  }

  Frequency lowestFrequency() { ret oscillators.firstKey(); }

  Frequency highestFrequency() { ret oscillators.lastKey(); }

  int minWindowSize() {

    ret iceil(sampleRate*lowestFrequency().interval());

  }

  S stats() { ret renderVars(+audio, +haarFeaturesChecked); }

  // how much horizontal area have we covered

  L<Double> coverageByFreq() {

    ret mapReversed(oscillators, o -> totalLengthOfDoubleRanges(keys(o.intensities)));

  }

  IIntegralImage imageColumn(DoubleRange timeRange) {

    int n = nFrequencies();

    double[] col = rawImageColumn(timeRange);

    double[] col2 = normalizeDoubles(col, 255);

    ret bwIntegralImageFromFunction(1, n, (x, y) -> ifloor(col2[y]));

  }

  double[] rawImageColumn(DoubleRange timeRange) {

    int n = nFrequencies();

    var l = reversedList(oscillators);

    double[] col = new[n];

    for y to n: {

      TreeMap<DoubleRange, Channels<Complex>> l2 = l.get(y).intensities;

      // find right period

      DoubleRange period = l2.floorKey(timeRange);

      if (period == null) continue;

      double amplitude = l2.get(period).first().abs();

      col[y] = amplitude;

    }

    ret col;

  }

  MakesBufferedImage image(DoubleRange r, double pixelsPerSecond) {

    new L<IIntegralImage> images;

    double len = r.length()/sampleRate;

    int pixels = iround(len*pixelsPerSecond);

    for i to pixels: {

      DoubleRange range = doubleRange(

        i/pixelsPerSecond*sampleRate+r.start(),

        (i+1)/pixelsPerSecond*sampleRate+r.start());

      //print(+range);

      images.add(imageColumn(range));

    }

    ret mergeBWImagesHorizontally(map toBWImage(images), spacing := 0);

  }

  double[] toAudio(IntRange sampleRange, Iterable<Oscillator> oscillatorsToRender default oscillators) {

    double[] samples = new[l(sampleRange)];

    class PerFreq {

      Oscillator o;

      DoubleRange period;

      Channels<Complex> intensity;

      *(Oscillator *o) {

        setPeriod(o.intensities.floorKey(toDoubleRange(sampleRange)));

        //printVars("PerFreq", f := o.frequency(), +sampleRange, +period, firstPeriod := firstKey(o.intensities));

      }

      void setPeriod(DoubleRange period) {

        this.period = period;

        intensity = period == null ?: o.intensities.get(period);

      }

      void nextPeriod(double t) {

        while (period != null && t >= period.end)

          setPeriod(o.intensities.higherKey(period));

      }

    }

    L<PerFreq> perFreqs = map(oscillatorsToRender, o -> new PerFreq(o));

    for i over samples: {

      double t = sampleRange.start+i; // t is time in samples

      double sum = 0;

      for (pf : perFreqs) {

        pf.nextPeriod(t);

        bool shouldPrint = toAudio_shouldPrint(i, l(samples));

        if (pf.intensity != null) {

          sum += renderer().renderPeriodSample(pf.period, pf.o.frequency(), pf.intensity.first(), t, shouldPrint);

        } else

          if (shouldPrint)

            printVars(f := pf.o.frequency(), +i, interval := pf.o.interval, +t, c := 0);

      }

      samples[i] = sum;

    }

    ret samples;

  }

  RenderFrequencySample renderer() {

    ret renderUsingSine

      ? new RenderFrequencySample_Sine

      : new RenderFrequencySample_SquareWave;

  }

  swappable bool toAudio_shouldPrint(int i, int n) { false; }

  void subtractedAudioSanityCheck() {

    if (rendered == null) ret;

    int n = l(rendered);

    for i to n: {

      reimmutable i;

      assertDoublesVeryClose_warnOnly(-> "Sample " + i,

        subtractedAudio.sampleSum(0, i, i+1),

        audio.sampleSum(0, i, i+1)-rendered[i],

        1e-5);

    }

  }

}