sclass ContinuousOscillators_TestAudioRestoration {
  S audioURL = "https://mouth.gazelle.rocks/beaHTML/15/5502%2Frecognize+-+61351880347bd84b28c7f42d779cde94.mp3";

  short[] inputSamples;
  double inputSampleRate = 48000;

  settable bool separateFrequencies;
  settable bool cumulative;
  settable bool separateRealAndImag;
  settable bool playRestored;
  settable bool showImage;
  bool useInputSampleRate = true;
  double silence = 0.2;
  double minFreq = 0, maxFreq = 48000;
  int useEveryNthOscillator = 1;
  
  short[] resampled;
  new L<double[]> audioOutput;
  File wavOut;
  
  new ContinuousOscillators co;

  void loadAudioURL {
    File mp3 = loadBinaryPageOnce(audioURL);
    loadMP3(mp3);
  }
  
  void loadMP3(File mp3) {
    printFileInfo(+mp3);
    File wav = mp3ToWAVUnlessExists(mp3, replaceExtension(mp3, ".mp3", ".wav"));
    loadWAV(wav);
  }
   
  void loadWAV(File wav) {
    //playWAV(wav);
    printFileInfo(+wav);
    wavOut = appendToBaseFileName(wav, ".restored");
    
    temp WAVDecoder decoder = new(wav);
    print(sampleRate := decoder.sampleRate);
    int n = decoder.totalValues();
    inputSamples = decodeWAVToMonoSamples(decoder, n);
    inputSampleRate = decoder.sampleRate;
  }
  
  // assumes soundSource is using the default inputSampleRate
  void setAudio(double seconds, VF1<double[]> soundSource) {
    inputSamples = concatShortArrays(soundSourceToShortArrays(seconds, soundSource, 1));
    print(inputSamples := sfu(takeFirst(10, inputSamples));
  }
  
  run {
    if (useInputSampleRate) co.sampleRate = inputSampleRate;
    resampled = convertSampleRate_shortArray_simple(inputSampleRate, co.sampleRate, inputSamples);
    print(samples := l(resampled));
    
    IQuerySound sound = shortSamplesToQuerySound(resampled);
    
    //co.setAudio(new StaticAudioSample(ll(resampled), 1, co.sampleRate));
    co.makeFullWindow(sound, intRange(0, l(resampled)));
    print(windowBounds := co.windowBounds());
    co.startOscillators();
    co.stepTo(co.windowBounds());
    
    print(coverageByFreq := co.coverageByFreq());
    
    restoreAudio();
    double[] audio = concatDoubleArrays(audioOutput);
    double[] audioNorm = normalizeDoubles(audio, 32767);
    
    if (wavOut != null) {
      shortArrayToMonoWAVE(doubleToShortArray_iround(audioNorm), wavOut, co.sampleRate);
      printFileInfo(+wavOut);
      if (playRestored) playWAV(wavOut);
    }
    
    print(audioNorm := sfu(takeFirst(10, audioNorm));
    
    double pixelsPerSecond = 100;
    if (showImage)
      dm_showPixelatedImage(co.imageForWindow(pixelsPerSecond));
  }
  
  void restoreAudio {
    L<ContinuousOscillators.Oscillator> oscillators = new L;
    //double ratio = 2; // octave
    int nOscillators = l(co.oscillators);
    var lOsc = asList(co.oscillators);
    
    for i to nOscillators: {
      if ((i % useEveryNthOscillator) == 0)
        oscillators.add(lOsc.get(i));
    }
    
    oscillators = filter(oscillators, o -> between(o.frequency()!, minFreq, maxFreq));
    
    print("Using " + nFrequencies(oscillators) + ": " + map(oscillators, o -> formatDouble(o.frequency()!, 1));
    
    if (separateFrequencies) {
      wavOut = appendToBaseFileName(wavOut, ".separate");
      double[] cum = null;
      for (o : oscillators) {
        pnl(map(o.intensities, (k, v) ->
          o.frequency() + " " + formatDouble(k, 1) + " = " + renderComplexWithAngle(div(v.first(), l(k)*32767))));
        co.toAudio_shouldPrint = (i, n) -> i < 30 || i >= n-10;
        for (xx : separateRealAndImag ? ll(1, 2) : ll(0)) {
          co.realOnly = (xx & 1) != 0;
          co.imagOnly = (xx & 2) != 0;
          double[] newAudio = co.toAudio(ifloor(co.windowBounds()), ll(o));
          if (cumulative) cum = newAudio = doubleAdd(cum, newAudio);
          audioOutput.add(newAudio); 
          audioOutput.add(repDouble(iround(co.sampleRate*silence), 0));
          co.toAudio_shouldPrint = null;
          co.realOnly = co.imagOnly = false;
        }
      }
    }
      
    audioOutput.add(co.toAudio(ifloor(co.windowBounds()), oscillators));
  }
}