Signale vergleichen

[Horsti (Gast)]

Hallo,
ich habe mit verschiedenen Sensoren eine Sinussignal gemessen. Dieses 
hat immer die gleiche Frequenz, gleiche Länge und gleiche Abtastrate. 
Jetzt möchte ich auf basis von dem genausten Sensor, die Anderen in 
Bezug auf ihrer Genauigkeit vergleichen. Ich habe schon versucht jeweils 
die Kreuzkorrelation mit Matlab zu berechnen, aber bin nicht in der Lage 
das Resultat zu deuten. Woran erkenn ich das die Messung verrauscht oder 
verzerrt ist? Wie würdet ihr die relative Genauigkeit bestimmen?

Hab mal ein Beispiel von der Kreuzkorrelation angehängt.

[---- (Gast)]

Diese Informationen können direkt aus dem Spektrum heraus bestimmt 
werden. Dazu gibt es freie Tools wie z.B. VisualAnalog von Analog 
Devices (Fensterung nicht vergessen!).
Gruss

[Matlab (Gast)]

Das was dir evtl. besser hilft ist der Korrelationskoeffizient:
http://de.wikipedia.org/wiki/Korrelationskoeffizient

Den kann auch Matlab berechnen.

[Matlab (Gast)]

Matlab schrieb:
> Den kann auch Matlab berechnen.

http://www.mathworks.de/help/techdoc/ref/corrcoef.html
http://www2.hs-esslingen.de/~mohr/mathematik/statistik/regression_mathe3_it.pdf
http://books.google.de/books?id=-nl5U_ZgoJUC&pg=PA149&lpg=PA149&dq=Korrelationskoeffizient+matlab&source=bl&ots=RREXZdJ4h7&sig=QMCm9Qn_WWnOtEj92BWJ82KlaS4&hl=de&sa=X&ei=Mhi-T4nSMMfY8QOpmKAh&ved=0CHIQ6AEwBQ#v=onepage&q=Korrelationskoeffizient%20matlab&f=false

[---- (Gast)]

Der Korrelationskoeffizient wird nicht viel bringen, denn er will wissen 
ob das Signal verrauscht oder verzerrt ist. Solche Grössen bestimmt man 
besten über den SNR und THD. Aus diesen Werten lässt sich direkt den 
genausten Sensor bestimmen.

[Horsti (Gast)]

Danke für die Antworten. Werde mir THD und SNR mal genauer anschauen. 
Dafür brauche ich ja das Frequenzspektrum. Mein Signal ist 500 Samples 
lang und wurde mit 7kHz abgetastet. Habe ich da nicht nur sehr 
beschränkte Möglichkeiten über eine FFT an das Spektrum zu kommen? Wenn 
ich ein Sinus mit 50 Hz analysiere bekomme ich den Peak im 
Frequenzbereich bei 52 Hz und die Harmonischen Oberschwingung sind nicht 
sichtbar. :(

[---- (Gast)]

Wen du die Y Achse in deinem Spektrum logarithmisch darstellst, sollten 
die Harmonischen sichtbar werden.

[Detlef _. (detlef_a)]

Horsti, Du hast bei der FFT Methode die Fenstereffekte mit drin, die 
machen Dir die Spektrallinie breit und erschweren die Bestimmung des SNR 
(SNR hat übrigens mit Spektrum nichts zu tun, das ist das Verhältnis vom 
Signal zum Rauschen). SNR Bestimmung geht auch ohne FFT, script anbei. 
Für schlechtere SNRs als 25dB oder so geht das allerdings nicht so gut, 
da muß man was anderes nehmen.

hope that helps
Cheers
Detlef


clear
w=15.3*2*pi/128
sig = 1.34*sin(w*(0:127)+0.123);
Psig=sum(sig.^2);
noise=1e-2*randn(size(sig));
Pnoise=sum(noise.^2);
sig=sig+noise;
Psn=sum(sig.^2);
n=length(sig);
%sig(n-1)+sig(n+1)=2*cos(w)*sig(n);
s1=sig(1:n-2)+sig(3:n);
s2=sig(2:n-1);
% Regression s1 auf s2
M= [s2.' ones(length(s1),1)];
coff=inv(M'*M)*M'*s1.';
wn=acos(coff(1)/2);

M=[cos((0:n-1)*wn).' sin((0:n-1)*wn).'];
coff=inv(M'*M)*M'*sig.';

cc=coff(2)+j*coff(1);
amp=abs(cc);
ph=angle(cc);

sigr=amp*sin(wn*(0:127)+ph);
noiser=sig-sigr;
Pnoiser=sum(noiser.^2);
Psigr=sum(sigr.^2);
disp(sprintf('SNR vorgegeben: %.2fdB SNR bestimmt  %.2fdB', ...
    10*log10(Psig/Pnoise),10*log10(Psigr/Pnoiser)));
n=length(sig);
plot([sig ; sigr].')

[---- (Gast)]

SNR lässt sich mit einem FFT sehr gut bestimmen, da sich die Signal und 
Rausch Leistungen einfach daraus ermitteln lassen.

Ich empfehle dir das VisualAnalog, da kannst du ausser der Fensterung 
nicht viel falsch machen.

[Horsti (Gast)]

@Detlef:
Danke, dass du dir die Mühe gemacht hast ein Beispielskript zu posten. 
Ehrlich gesagt verstehe ich nicht, was du für ein Verfahren benutzt. Hab 
dein Besipiel mal mit einem Signal ausprobiert, was dem von mir 
untersuchten ähnlich ist. Dabei ist die Abweicheung des bestimmten SNR 
zum vorgegeben zu groß. Kann natürlich sein, dass ich irgendeinen Mist 
gebaut habe. Wie gesagt bin nicht ganz durch gestiegen. Verändertes 
Skript ist im Anhang

@Anderen:
Bezüglich THD

Hat natürlich geholfen die Y Achse logarithmisch darzustellen. Dennoch 
bereitet mir der Leakage Effekt Probleme. Die harmonischen sind dadurch 
nicht sichtbar. Um das zuvermeiden wird empfohlen, die Anzahl der 
Samples so zuwählen, dass sie ein ganzahliges vielfaches der 
Periodendauer des Signals sind. Damit die Harmonischen genau in den 
Nullstellen des sin(k)/k sitzen.
Kein Problem.
f = 50 Hz, T= 0.02 s, Fs = 7 kHz, N = Anzahl der Samples, B = 
Beobachtungsfenster
B = n * T= 0.06 = N * 1/7kHz
n = 3, N = 42
Das Problem die Nullstellen liegen eben nicht genau auf den Vielfachen 
der Signalfrequenz

Ok zweite Möglichkeit wäre den Leakageffekt mit einer Fensterung zu 
verringern. Aber wie geht das? Nehm ich das Fenster und verschieb es 
immer um ein ganzahliges der Signalfrequenz um den Betrag der 
verschiedenen Harmonischen zu bekommen?

Visual Analog habe ich mir mal angeschaut. Komme damit aber überhaupt 
nicht klar. Bekomm noch nicht mal die csv Datei geladen. Außerdem würde 
ich gerne Matlab benutzen, damit die Grafiken zu Illustration konsistent 
sind.

Eine Weitere Frage zur Berechnung des SNR. Kann ich mit einem Tiefpass 
Filter die Leistung meines Signals und mit einem Hochpassfilter die 
Leistung des Rauschens bestimmen? Oder wie geht das?

[Werner (Gast)]

Horsti schrieb:
> Kann ich mit einem Tiefpass Filter die Leistung meines Signals und mit
> einem Hochpassfilter die Leistung des Rauschens bestimmen?

Nein, ein Rauschspektrum enthält i.A. alle Frequenzen. Mit einem Hochpaß
würdest du den niederfrequenten Anteil abschneiden.
In der Fouriertransformierten kannst du beides direkt ablesen.

Horsti schrieb:
> Wie würdet ihr die relative Genauigkeit bestimmen?

Was meinst du mit Genauigkeit? Verwechselst du da was mit 
Signal-Rausch-Verhältnis?

[Horsti (Gast)]

Wie bekomme ich denn dann die Rauschleistung?
Stimmt ich bringe, dass immer durcheinander. Ich möchte bestimmen, 
welcher Sensor am besten ist. Also hohes SNR, niedriges THD und nochwas?

[---- (Gast)]

Habe vor Jahren dieses Skript für die Spezifikation eines ADCs 
geschrieben, daraus sollte ersichtlich sein wie so was gemacht wird.

% Dieses Script analysiert die Signaleigenschaften wie, snr, thd, sinad, sfdr 
% und enob.
load JitterCorrectedData_0.txt;
signal = JitterCorrectedData_0;
numberOfSamples = length(signal);
samplingRate = 5e9;
maxNumberOfHarmonics = 5;
harmonicWidth = 9; % breite der Peaks

% Spektrum berechen
window = blackman(numberOfSamples);
signal = signal(:,1);
signal = signal .* window;
signal = abs(fft(signal));
signal = signal(1:numberOfSamples/2);

% DC und erste Harmonische suchen. Der grösste Wert, ausser DC, wird als erste 
% Harmonische gedeutet.
levels = signal(1);
positions = 0;
[level positionH1] = max(signal(2:numberOfSamples/2));
positionH1 = positionH1;
levels = [levels level];
positions = [positions positionH1];

% Alle anderen Harmonischen finden.
currentPosition = 2 * positionH1;
while currentPosition-2 < numberOfSamples/2 && length(levels) < maxNumberOfHarmonics
  [level lastPosition] = max(signal(currentPosition-2:min(currentPosition+2, numberOfSamples/2)));
  lastPosition = lastPosition + currentPosition - 3;
  levels = [levels level];
  positions = [positions lastPosition];
  currentPosition = currentPosition + positionH1;
end;

% Pegels skalieren damit die 1. Harmonische eins beträgt.
levelScalingFactor = 1 / levels(2);
levels = levels * levelScalingFactor;
signal = signal * levelScalingFactor;

% Noise vom Signal trennen und Leistungen der Harmonischen berechnen.
noiseMasked = signal;
harmonicWidth2 = ceil(harmonicWidth/2);
noiseMasked(1:harmonicWidth2) = zeros(1,harmonicWidth2)*-1;
powers = sum(signal(1:harmonicWidth2).^2);
for i=2:length(positions)
  from = positions(i) - harmonicWidth2 + 2;
  to = from + harmonicWidth - 1;
  noiseMasked(from:to) = zeros(1,harmonicWidth)*-1;
  power = sum(signal(from:to).^2);
  powers = [powers power];
end;

% Lehrstellen aus dem noiseMasked Vektor entfernen.
noiseContinuous = [];
for i=1:length(noiseMasked)
  if noiseMasked(i) ~= 0
    noiseContinuous = [noiseContinuous; noiseMasked(i)];
  end;
end;

% snr, thd, sinad, sfdr und enob berechnen.
noisePower = sum(noiseContinuous.^2);
signalPower = powers(2);
harmonicsPower = sum(powers(3:end));

snr = 10*log10(signalPower/noisePower);
thd = 10*log10(harmonicsPower/signalPower);
sinad = 10*log10(signalPower/(noisePower+harmonicsPower));

[highestHarmonicsLevel highestHarmonicsPosition] = max(levels(3:end));
[highestNoiseLevel highestNoisePosition] = max(noiseMasked);
sfdr = mag2db(levels(2)/max([highestHarmonicsLevel highestNoiseLevel]));
enob = (sinad - 1.76) / 6.02;

% Resultate ausgeben.
disp(['SNR    = ' num2str(snr) ' dB']);
disp(['THD    = ' num2str(thd) ' dB']);
disp(['SINAD  = ' num2str(sinad) ' dB']);
disp(['SFDR   = ' num2str(sfdr) ' dB']);
disp(['ENOB   = ' num2str(enob) ' bit']);

% Plotten.
xScalingFactor = samplingRate/numberOfSamples;
x = linspace(0,samplingRate/2-xScalingFactor,numberOfSamples/2);
plot(x, mag2db(signal));
hold on
plot(x, mag2db(noiseMasked), 'LineWidth', 3);
plot(0, mag2db(levels(1)), 'rs');
plot(positions(2)*xScalingFactor, mag2db(levels(2)), 'rv');
plot(positions(3:end)*xScalingFactor, mag2db(levels(3:end)), 'ro');
plot((highestNoisePosition-1)*xScalingFactor, mag2db(highestNoiseLevel), 'rd');
hold off;

legendText = num2str(round(mag2db(1/levels(3))));
for i=4:length(levels)
  legendText = [legendText  ', ' num2str(round(mag2db(1/levels(i))))];
end;

legend('signal', ...
       'noise', ...
       ['dc component ' num2str(round(mag2db(1/levels(1)))) ' dBc'], ...
       'fundamental signal', ...
       ['harmonics ' legendText ' dBc'], ...
       ['highest spur ' num2str(round(mag2db(1/highestNoiseLevel(1)))) ' dBc']);

title('Frequency Domain Analysis');
xlabel('Frequency [Hz]');
ylabel('Level [dB]');
set(gca,'YLim',[-80 20]);

[Horsti (Gast)]

Vielen Dank. Mittlerweile ist mir die Berechnung von THD, SNR auch klar 
geworden. Aber vielleicht werde ich auch SFDR und ENOB mit in die 
Analyse mit einbeziehen. Obwohl sich ENOB, wenn ich das richtig 
verstanden habe, mehr auf ADCs bezieht. Gibt es eine Fenstermethode die 
generell zu empfehlen ist? Mit hamming bekomme ich schon ganz gute 
Werte. Mit kaiser oder blackmanharris bekomme ich qualitativ die selben 
Ergebnisse. Dafür ist das SNR generell etwas niedriger.