[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

[sc-users] Spectral subtraction



Hi there folks!
A while ago I asked if anyone had any tip on subtracting one spectra from
another. This is what I came up with. ( would be nice if you could select
function for the input and oitput windows in the IFFT.) ex.
IFFT.ar(fftsize, fftoffsets, cosineTable, filterWindow, nil, fft.real,
fft.imag, inputWindowFunction: [mul, add, sub, div].choose);
to be able to subtract add, or divide spectras.
I wanted to subtract a fofvoice singing ah from a real voice singing the
same ah to get the difference ( the humanity left when the static fofvoice
was subtracted) Of cource you could use it to ex. remove the spectra of a
symphony orchestra from a white noise or something. I searched lots of
applications but none seemed to to what I wanted. MetaSynth has a spectral
filter but multiplying. If anyone has comments for improvement I'd be glad
to take part of it; I guess that the diff of the magnitudeDiff and the
result is due to convolution with the window?

Jesper


(  // 30/12-2002 commented version for SC list
      var fftsize, numffts, fftoffsets, window, cosineTable, input,
inputAmp, threshhold, gate, filterWindow;
    var src1, src2, fft1, fft2, fft3, ifft, out, max, ind, voiceFilename,
fofFilename, voiceSoundFile, fofSoundFile, voiceSignal, fofSignal, size;
    var voiceImag, voiceTable, voiceReal, voiceFft, voiceMagnitude,
plotsize, fofImag, fofTable, fofReal, fofFft, fofMagnitude, magnitudeDiff,
freq, hasFreq;
    var magnitudeDiffPlot, filterWindowPlot, voiceMagnitudePlot,
fofMagnitudePlot, w, magnitudeDiffTable, filterWindowTable,
voiceMagnitudeTable, fofMagnitudeTable;
    var filterWindow2, result, resultPlot, resultPlotTable, min, s1, s2, s3,
s4, s5, maxFreq, zoom;
    var resultMagnitude, resultFft, resultImag, resultTable, resultReal,
resultMagnitudePlot, resultWindowTable;
    var fftWindow, diffPeak, resultPeak, peak;
    
    w = GUIWindow.new("Spectral Subtraction", Rect.newBy(5, 41, 1016,
724)).backColor_(rgb(120,120,120));
    voiceMagnitudeTable = WavetableView.new( w, Rect.newBy(10, 25, 999,
110), nil, -0.1, 1.1);
      fofMagnitudeTable = WavetableView.new( w, Rect.newBy(10, 165, 999,
110), nil, -0.1, 1.1);
     magnitudeDiffTable = WavetableView.new( w, Rect.newBy(10, 305, 999,
110), nil, -0.1, 1.1);
      filterWindowTable = WavetableView.new( w, Rect.newBy(10, 445, 999,
110), nil, -0.1, 1.1);
      resultWindowTable = WavetableView.new( w, Rect.newBy(10, 585, 999,
110), nil, -0.1, 1.1);
    
    s1 = StringView.new( w, Rect.newBy(112, 0, 800, 20),   "
Real Voice FFT Magnitude:").labelColor_(Color.white);
    s2 = StringView.new( w, Rect.newBy(112, 138, 800, 20), "
Fof Voice FFT Magnitude:").labelColor_(Color.white);
    s3 = StringView.new( w, Rect.newBy(112, 278, 800, 20), "       Magnitude
Difference (Real Voice FFT Magnitude - Fof Voice FFT Magnitude):
").labelColor_(rgb(200,200,80));
    s4 = StringView.new( w, Rect.newBy(112, 418, 800, 20), "Filter
(Magnitude Difference / Real Voice Spectra) to Muliply with Real Voice
Spectra:").labelColor_(rgb(200,200,80));
    s5 = StringView.new( w, Rect.newBy(112, 558, 800, 20), "
Resulting Spectra of multiplying Real Voice with
Filter:").labelColor_(rgb(255,255,0));
    
    numffts = 4;   
    fftsize = 4096;
    
    
    // Uncomment and Exchange for your own soundfiles
    /*
    voiceFilename = "Macintosh HD:Program:SC2.2.15 f:woodwinds:Fem Voices G3
L 1 sec.aiff"; 
    fofFilename = "Macintosh HD:Program:SC2.2.15 f:woodwinds:snd5.aiff";
        
    voiceSoundFile = SoundFile.new;
    if (voiceSoundFile.read(voiceFilename), { voiceSignal =
voiceSoundFile.data.at(0); });
    
    fofSoundFile = SoundFile.new;
    if (fofSoundFile.read(fofFilename), { fofSignal =
fofSoundFile.data.at(0); });
    */
    
    
    // test signals instead of sound files
    voiceSignal = Synth.collect({ Blip.ar( 86.1328, 24) }, fftsize /
Synth.sampleRate); 
    fofSignal = Synth.collect({ Blip.ar( 86.1328, 20) }, fftsize /
Synth.sampleRate); 
    
    fftWindow = Signal.hanningWindow( fftsize ); // window to reduce
artifacts 
    
    // use only first 4096 samples of Voice and get a normalized magnitude
signal    
    voiceReal = voiceSignal.copyRange(0, fftsize - 1) * fftWindow;
    voiceImag = Signal.newClear(fftsize);
    voiceTable = Signal.fftCosTable(fftsize);
    voiceFft = fft(voiceReal, voiceImag, voiceTable);
    voiceMagnitude = voiceFft.magnitudeApx.normalize;
    
    // use only first 4096 samples of fofVoice and get a normalized
magnitude signal   
    fofReal = fofSignal.copyRange(0, fftsize - 1) * fftWindow;
    fofImag = Signal.newClear(fftsize);
    fofTable = Signal.fftCosTable(fftsize);
    fofFft = fft(fofReal, fofImag, fofTable);
    fofMagnitude = fofFft.magnitudeApx.normalize;
    
    magnitudeDiff = (voiceMagnitude.copy - fofMagnitude.copy).thresh(0.0);
// ideal result
    filterWindow = (magnitudeDiff.copy / voiceMagnitude.copy).thresh(0.0);
// calculate filter
    // the result of multyplying filterWindow with Voice should result in
magnitudeDiff
    // If not, should be convolution with window
    
    diffPeak = magnitudeDiff.peak; // get peak of magnitudeDiff
    
    plotsize = (fftsize / 16).asInteger - 1; // calculate how many samples
to plot; 
    maxFreq = ((1 + plotsize) * (Synth.sampleRate / fftsize)); // calculate
maxFreq of plot
    
    s1.label = s1.label ++ " First " ++ (1 + plotsize).asString ++ "
Frequencys " ++ " ( 0 to " ++ maxFreq.asString ++ " htz )";
    s2.label = s2.label ++ " First " ++ (1 + plotsize).asString ++ "
Frequencys " ++ " ( 0 to " ++ maxFreq.asString ++ " htz )";
    s3.label = s3.label ++ " First " ++ (1 + plotsize).asString ++ "
Frequencys " ++ " ( 0 to " ++ maxFreq.asString ++ " htz )";
    s4.label = s4.label ++ " First " ++ (1 + plotsize).asString ++ "
Frequencys " ++ " ( 0 to " ++ maxFreq.asString ++ " htz )";
    s5.label = s5.label ++ " First " ++ (1 + plotsize).asString ++ "
Frequencys " ++ " ( 0 to " ++ maxFreq.asString ++ " htz )";
    
    // get first plotsize samples of magnitudes
    magnitudeDiffPlot = magnitudeDiff.copyRange(0, plotsize); // parts to
plot
    voiceMagnitudePlot = voiceMagnitude.copyRange(0, plotsize);
    fofMagnitudePlot = fofMagnitude.copyRange(0, plotsize);
    filterWindowPlot = filterWindow.copyRange(0, plotsize);
    // Convert to wavetables and plot
    magnitudeDiffTable.wavetable = magnitudeDiffPlot.asWavetable; // change
to wavetables
    voiceMagnitudeTable.wavetable = voiceMagnitudePlot.asWavetable;
    fofMagnitudeTable.wavetable = fofMagnitudePlot.asWavetable;
    filterWindowTable.wavetable = filterWindowPlot.asWavetable;
    // for playback
    fftoffsets = Array.series(numffts, 0, fftsize/numffts);
    window = Signal.hanningWindow(fftsize, fftsize/2);
    cosineTable = Signal.fftCosTable(fftsize);

result = // snapshot of result
Synth.collect({ 
    // src1 = Warp1.ar( voiceSignal, 0.0 );
    src1 = PlayBuf.ar( voiceSignal, Synth.sampleRate,  1, 0, 0, 44100);
    fft1 =  FFT.a2sr(fftsize, fftoffsets, cosineTable, window, nil, src1,
0.0); 
    ifft = IFFT.s2ar(cosineTable, filterWindow, nil, fft1.real, fft1.imag);
    out = Mix.ar(ifft.real);
    out }, 
    ( fftsize * 2) / Synth.sampleRate
    );
    
    // compencating for fft delay
    resultReal = result.copyRange(fftsize, fftsize + fftsize - 1) *
fftWindow;  
    resultImag = Signal.newClear(fftsize);
    resultTable = Signal.fftCosTable(fftsize);
    resultFft = fft(resultReal, resultImag, resultTable);
    // get normalized magnitude of result signal
    resultMagnitude = resultFft.magnitudeApx.copy.normalize;
    resultPeak = resultMagnitude.peak; // get peak of result magnitude (
should be 1 of cource )
    // get difference of peaks
    peak = diffPeak / resultPeak;
    // get first plotsize samples of magnitudes
    resultMagnitudePlot = resultMagnitude.copyRange(0, plotsize) * peak;
    // Convert to wavetables and plot
    resultWindowTable.wavetable = resultMagnitudePlot.asWavetable;
    
    // play the result of subtracting fofVoice from real Voice
Synth.play({ arg synth;
    src1 = PlayBuf.ar( voiceSignal, Synth.sampleRate,  1, 0, 0,
voiceSignal.size);
    // If you have Chad's excellent warp class you can play a nonstatic loop
of the signal 
    // src1 = Warp1.ar( voiceSignal, 0.0 );
    fft1 =  FFT.ar(fftsize, fftoffsets, cosineTable, window, nil, src1,
0.0); 
    ifft = IFFT.ar(fftsize, fftoffsets, cosineTable, filterWindow, nil,
fft1.real, fft1.imag);
    out = Mix.ar(ifft.real);
    // add some reverb
    4.do({ out = AllpassN.ar( out, 0.05, [ 0.01 + 0.05.rand, 0.01 +
0.05.rand ], 1 ) });
    out
    });
  // w.close;    
)