Wieso wird die Kontaktierung bei 2,98s nicht erkannt und die beiden benachbarten schon?
Hier ist der Auschnitt. Man erkennt das vor dem eigentlichen Event etwas stärkeres Rauschen ist, wo die hellblaue Kurve schon ansteigt, das ist vor den anderen Events nicht so gewesen. Das müsste mit einem breiteren Fenster (stärkere Filterung) besser gehen. *ausprobier*
Sieht doch schon besser aus.
Ja, du hast einiges an dem Auswerteprogramm verändert.
Aber ohne sich in die Programmiersprache einzulesen, ist es kaum nachvollziehbar.
Ok ich versuche den Code mal zu erklären.
Erstmal allgemeine Syntaxgrundlagen.
Mit einem Semikolen hinter einem Befehl kann der Output unterdrückt werden und ein Prozentzeichen leitet einen Kommentar ein. Hilfe zu einem Befehl bekommt man schnell mit help <Befehlname>.
Code:
close all; % alle Plot-Fenster schließen
[data,fs]=wavread('Reed-Kontakt mit DMM.wav'); % in data steckt das Signal, fs ist die Abtastrate
U=1.5; % Reifenumfang in Metern
s=6.5; % Schwellwert für die Erkennung der Events
Um die Peaks zu erkennen, muss man den Grundpegel kennen. Dafür verwende ich einen laufenden Mittelwert. Das Fenster darf dabei nicht zu breit sein, damit es nicht mehrere Peaks beinhaltet. Es sollte aber auch nicht zu kurz sein, da Rauschen sonst eine größere Rolle spielt. Ursprünglich hatte ich die Länge eines Events als orientierung genommen, damit ich auch zwei Events, die kurz hintereinander kommen, erkennen kann. Da das Event selbst ein abklingenden Verlauf hat, sollte man sogar zwei direkt aufeinander folgende Events erkennen können. Das ist allerdings fern jeder Praxis, da wir so hohe Geschwindigkeiten nicht messen werden. Also habe ich mal eine maximale Geschwindigkeit von 200km/h angesetzt (ich habs auch mal auf 60 gesetzt, das war aber zu viel des guten und 200 resultierte in guten Ergebnissen ...
)
Code:
vmax=200; % theoretische maximale erkennbare Geschwindigkeit
eventlength=U*fs/vmax*3.6;
ww = eventlength/2; %Fensterbreite
Bei 200km/h und der aktuellen Abtastrate würde man alle U*fs/vmax*3.6 Datenpunkte einen Peak erwarten. Für das Fenster nehme ich daher die halbe Anzahl an Datenpunkten.
Code:
if mod(ww/2,2) ~= 0
ww = ww -1;
end
Das Fenster soll eine gerade Anzahl an Datenpunkten haben. Das ist im weiteren Verlauf wichtig.
Code:
absdata=abs(data); % Betrag des Signals
w=[zeros(1,ww) ones(1,ww)./ww]; % Fenster
c=conv(absdata,w); % Berechnung des laufenden Mittelwerts
conv steht für
convolution und bietet mir die Möglichkeit den laufenden Mittelwert schnell zu berechnen. Dabei wird das Fenster quasi über das Signal geschoben und die Summe der Produkte gebildet. Wir interessieren uns nur für die Vergangenheit, deswegen besteht das Fenster zur Hälfte aus Nullen und wir wollen den Mittelwert haben, deswegen ist die Summe der einzelnen Fenstereinträge gleich eins.
Code:
t=[1:length(absdata)]./fs; % Zeitvektor generieren
c=c(floor(0.5*length(w)):end-floor(0.5*length(w)));
Durch die Berechnung mit conv wird der Ergebnisvektor länger als das ausgangssignal (Das Fenster wird über den Anfang und das Ende Hinaus geschoben.) Diese zusätlichen Werte interessieren uns nicht und werden mit dem c=c(... hier entsorgt.
Code:
tevt=t(absdata./c >6.5); % Eventerkennung
Ok hier steckt tatsächlich sehr viel drin
t ist der Zeitvektor und mit t(1) kann ich mir den ersten Wert ansehen (mit t(end) den letzten und mit t(20:30) die 11 Werte 20,21,...,30).
./ ist ein Elementoperator, dafür müssen die beiden Argumente (hier absdata und c) die gleiche Länge haben. Es wird hier also jedes Element von absdata durch das entsprechende Element von c geteilt. Also t(1)/c(1) und t(10)/c(10) und so weiter.
Wenn das Signal in absdata nun größer ist als das 6.5 fache von c wird der dazugehörige Zeitpunkt aus t in tevt geschrieben.
Code:
tevt([0 diff(tevt)<ww/fs]>0)=[];
Da da ein Event nicht nur aus einer halben Sinuswelle besteht sondern aus vielen abklingenden. Wird zwangsläufig hier und da ein Event quasi doppelt erkannt. Da wir aber wissen, dass innerhalb eines Fensters nur ein Event sein kann (So haben wir das Fenster schließlich entworfen.), können wir die falschen Events erkennen.
diff(tevt) Berechnet die Abstände der erkannten Events. Wenn diese Abstände kleiner als die Fensterbreite sind löschen wir den späteren Wert.
diff(tevt) gibt also einen Vektor zurück der ein Element weniger hat als tevt. durch [0 diff(tevt)] erhalten wir einen Vektor mit der gleichen Länge. Der Vektor "diff(tevt)<ww/fs" besteht nur aus Nullen und Einsen die für Falsch und Wahr entsprechend der kleiner als Bedingung stehen. Durch das Erweitern des Vektors mit der Null am Anfang wird immer der zweite Wert gelöscht, denn wir Wollen ja den Anfang des Events erkennen.
Code:
figure();
plot(t,absdata,'b-',t,c,'c-',tevt,0,'r*')
figure();
v=U./diff(tevt)*3.6;
tv=zeros(length(v),1);
for i=1:length(tevt)-1
tv(i) = mean(tevt(i:i+1));
end
plot(tv,v);
xlabel('Zeit [s]');ylabel('Geschwindigkeit [km/h]');
Wenn ich die Geschwindigkeit durch zwei Zeitpunkte und den in dieser Zeit zurückgelegten Weg berechne, ist das eine Durchschnitsgeschwindigkeit und der dazugehörige Zeitpunkt liegt genau zwischen den beiden Zeitpunkten. Das mache ich mit der for-Schleife. Sonst wird hier der Kram nur dargestellt.
Man kann jeden Befehl einfach so in Oktave eingeben oder das Skript mit dem Skriptnamen aufrufen, wenn es sich um Suchpfad befindet.
Einfach einbisschen rumspielen und ausprobieren. Und bei Fragen weiß google eigentlich immer weiter. Ein bisschen weiter weiß Google, wenn man nicht nach octave sondern nach matlab und dem entsprechenden Befehl oder der gewünschten Funktionalität sucht und erst im Anschluss schaut, ob der Befehl auch bei octave funktioniert oder was für eine Alternative es gibt.
Viele Grüße,
Nepomuk
