Circular spectrum visualizer

[Delphi.Narium]

Oder kann es doch noch irgendwo zu 'nem Überlauf kommen.

Eigentlich hätte ich nach den bisherigen Analysen immer mit Werte rund um 0 gerechnet (so wie im zweiten Screenshot).

Die Werte im ersten Screenshot liegen irgendwie alle weit außerhalb des Bereiches, den man in 'nem Integer unterbringen kann.

Solange es sich nur um positive Werte handeln sollte, also > 0, würden sie so gerade eben in 'nen Cardinal passen.

Da es aber (wie mir scheint) auch negative Werte geben kann, muss es wohl Int64 werden.

[Michael II]

Hoi EWeiss

die von dir genutzte uSpectrum.FFT Funktion rechnet falsch.

Nimm doch eine hier aus dem Forum, zum Beispiel diese hier


Ich habe die Unit aus dem Forum etwas gekürzt (Code hier unten). Diese Unit fügst du zu deinem Projekt hinzu:

Delphi-Quellcode:

unit uDFT;

interface

uses Math;

//
// Autor Matze - siehe: https://www.delphipraxis.net/597828-post1.html
//


type
  TComplex = record
    re, im: Extended;
  end;

  TComplexArray = array of TComplex;

procedure DFT(var a: TComplexArray);

implementation

function AddC(a, b: TComplex): TComplex;
begin
  Result.re := a.re + b.re;
  Result.im := a.im + b.im;
end;


function SubC(a, b: TComplex): TComplex;
begin
  Result.re := a.re - b.re;
  Result.im := a.im - b.im;
end;


function MulC(a, b: TComplex): TComplex;
begin
  Result.re := a.re * b.re - a.im * b.im;
  Result.im := a.re * b.im + a.im * b.re;
end;

function MakeC(re, im: extended): TComplex;
begin
  Result.re := re;
  Result.im := im;
end;


procedure shuffle(var a: TComplexArray; n, lo: Integer);
var I, m: Integer;
    b: TComplexArray;
begin
  m := n shr 1;
  setlength(b, m);
  for I := 0 to m - 1 do
    b[i] := a[lo + i];
  for I := 0 to m - 1 do
    a[lo + i + i + 1] := a[lo + i + m];
  for I := 0 to m - 1 do
    a[lo + i + i] := b[i];
end;


procedure DoFFT(var a: TComplexArray; n, lo: Integer; w: TComplex);
var I, m: Integer;
    z, v, h: TComplex;
begin
  if n and (n - 1) = 0 then
  begin
    if n > 1 then
    begin
        m := n shr 1;
        z := MakeC(1, 0);
        for I := lo to lo + m - 1 do
        begin
            h := SubC(a[i], a[i + m]);
            a[i] := AddC(a[i], a[i + m]);
            a[i + m] := MulC(h, z);
            z:=MulC(z,w);
        end;
        v := MulC(w, w);
        DoFFT(a, m, lo, v);
        DoFFT(a, m, lo + m, v);
        shuffle(a, n, lo);
    end;
  end;
end;

procedure DFT(var a: TComplexArray);
begin
    DoFFT(a, length(a), 0, MakeC(cos(2 * Pi / length(a)),
        sin(2 * Pi / length(a))));
end;
end.

Wenn du die verwendeten Dateitypen (dein TComplex und Matzes TComplex) nicht anpassen magst, dann ersetze in uSpectrum.pas die FFT Funktion durch diese hier (Code unten). (uDFT unter uses hinzuzufügen.)

Delphi-Quellcode:

uses ….uDFT; 

…
…


procedure TSpectrum.FFT( var Dat : array of TComplex );
var a : uDFT.TComplexArray;
    i, n : integer;
begin
  n := length( Dat );
  setlength( a, n );

  for i := 0 to n-1 do
  begin
    a[i].re := Dat[i].r;
    a[i].im := Dat[i].i;
  end;

  DFT( a );

  for i := 0 to n-1 do
  begin
    Dat[i].r := a[i].re/n;
    Dat[i].i := a[i].im/n;
  end;
end;

Ich verwende für meine Programme eine iterative Version von Cooley und Tukey. Ich speichere dabei sämtliche Einheitswurzeln einmal in einer Tabelle ab und greife dann auf diese zu. Eine solche iterative Lösung ist bei der von dir gewählten Problemgrösse FFFTSize=2048 aber nur ca. 6 Mal schneller.

Selbst auf meinem langsamen Notebook benötigt obige Lösung nur ca. 2/1000 Sekunden.


Viel Spass beim Codieren.

[EWeiss]

@Michael Danke werde es mal testen.
Ich habe jetzt aber das Problem Gain wird in der FFT nicht berücksichtigt.

Delphi-Quellcode:

  if FView = 0 then
    sr := (4096 * FGain) / FFFTSize
  else
  sr := 1 / FFFTSize;

  for i := 0 to (FFFTSize div 2) - 1 do
  begin
    Dat[i].r := LimitedSingleValue(Dat[i].r * sr);
    Dat[i].i := LimitedSingleValue(Dat[i].i * sr);
  end;

Wo muss ich eine vergleichbare Berechnung zur FFT hinzufügen?
Mit deinen neuen FFT habe ich vielleicht 10 Pixel bei normaler Visualisierung.

Mit der von Matze wird die Visualisierung fast korrekt angezeigt aber mit einer Auslastung von 80% eines CPU Kerns.

gruss

[EWeiss]

Der FFT von Matze ist leider auch nicht korrekt aber destotrotz mein Ergebnis ist jetzt zu 98% richtig. Einen Tick zu wenig! CPU last 2%
Den Rest bekomme ich auch noch hin.
Bin also doch nicht ganz so blöd wie ich hier hingestellt werde brauche lediglich etwas mehr zeit!
Benötigte nur einen Funktionierenden FFT! Weil meine Mathe Kenntnisse dafür nicht ausreichen (niemand ist perfekt)

Ich habe jetzt Zeit investiert aber was zählt das schon die wird mit unter für weniger sinnvolles verplempert.

gruss

[Michael II]

Hoi EWeiss

schön klappt's.

DFT sollte funktionieren, ich habe die Funktion anhand einiger Beispiele durchgerechnet und bei weniger plausiblen Argumenten mit den Werten meiner eigenen Funktion verglichen.

Zum View=0 Fall. Ich hatte übersehen, dass der Autor für diesen Fall die Funktion FFT "missbraucht". Sowas gehört eigentlich (wie er es beim Fall View=1 getan hat) in die Prozedur, welche die Werte grafisch darstellt.

Wenn du den View=0 Fall dennoch (wie im Original) in die FFT() einbauen willst, dann so:

Delphi-Quellcode:

procedure TSpectrum.FFT( var Dat : array of TComplex );
var a : uDFT.TComplexArray;
    i, n : integer;
    sr : single;
begin
  n := length( Dat );
  setlength( a, n );

  for i := 0 to n-1 do
  begin
    a[i].re := Dat[i].r;
    a[i].im := Dat[i].i;
  end;

  DFT( a );

  if FView = 0 then
  sr := FGain/50
  else
  sr := 1/n;

  for i := 0 to n-1 do
  begin
    Dat[i].r := a[i].re * sr;
    Dat[i].i := a[i].im * sr;
  end;
end;

Und falls du View nicht in FFT einbaust, dann musst du für View=0 die Spektren mit dem Faktor scale (Code unten) skalieren:

Delphi-Quellcode:

    0:
      begin
        scale := FFFTSize/FGain*50;

        for o := 0 to FOctaveCount - 1 do
        begin
          fl := True;
          q2 := q2 + OctAreaSize;

          i2 := i1 * 2;
          while i1 < i2 do
          begin
            b := scale*Sqr(Spectrum[i1].r * Spectrum[i1].r + Spectrum[i1].i * Spectrum[i1].i);

[EWeiss]

Danke dir erstmal wie immer

Ich habe es so gemacht.. mit einer anderen FFT die funktioniert sorry.

Delphi-Quellcode:

procedure TSpectrum.FFT(var Dat: array of TComplex);
var
  tp: array of TComplex;
  i: integer;
  sr: Extended;

begin

  SetLength(tp, FFFTSize);

  for i := 0 to FFFTSize - 1 do
  begin
    tp[i].re := Dat[i].re;
    tp[i].im := Dat[i].im;
  end;

  ForwardFFT(Dat, tp, FFFTSize);

  if FView = 0 then
    sr := (4096 * FGain) / FFFTSize
  else
  sr := 1 / FFFTSize;

  for i := 0 to (FFFTSize div 2) - 1 do
  begin
    Dat[i].re := tp[i].re * sr;
    Dat[i].im := tp[i].im * sr;
  end;

end;

Zitat:

Und falls du View nicht in FFT einbaust, dann musst du für View=0 die Spektren mit dem Faktor scale (Code unten) skalieren:

View ist nötig weil ich hier den Effekt umschalte zwischen Ring und Sector. Ohne funktioniert Rings nicht.

Ich kann deine FFT gerne nochmal gegentesten mit deiner Änderung.

gruss

[Michael II]

Ist klar, dass du View benötigst.

Ich würde es entweder wie im Code 2 gezeigt nachbessern; oder du kannst die angepasste FFT verwenden.