Beispiele zum Wavelet Explorer

ADDITIVE GmbH
Abt. Engineering & Science
Rohrwiesenstraße 2
D-61381 Friedrichsdorf/Ts.
  Tel.:  +49 6172 5905 30
  Fax:  +49 6172 77613
  email: mathematica@additive-net.de
  www:  http://www.additive.net.de

Die Mathematica Applications Library

Mathematica hat sich in den vergangenen Jahren als unentbehrliches Werkzeug im Ingenieurwesen etabliert. Es bietet dem professionellen Anwender, durch seinen mächtigen und fortschrittlichen Funktionsumfang, eine unerreichte Flexibilität bei der Analyse, Modellierung und Visualisierung der während der täglichen Ingenieurarbeit anfallenden Daten. Im Sommer 1994 startete Wolfram Research die Veröffentlichung der Mathematica Anwendungs-Pakete, die dem Benutzer sofort einsetzbare Funktionen aus speziellen Anwendungsbereichen bieten. Dem Anwender wird durch diese Pakete der schrittweise Einstieg in die Thematik bis zum Expertenwissen leicht gemacht. Die Zahl der Zusatzpakete wurde seit dieser Zeit auf über 30 Zusatzpakete erweitert.

Ein Übersicht zu den Zusatzpaketen ist aktuell im INTERNET zu finden.

http://www.additive-net.de/software/mathematica/packages/mma.packages.shtml

System Anpassung Demonstration

Mit den Beispielen in desem Dokument können Sie nur arbeiten, wenn Sie das Zusatzpaket installiert haben. Das ist bei Ihnen nicht der Fall.
Das Dokument stellt damit lediglich eine Demonstration der Möglichkeiten dar.
Nachfolgend werden Ihnen die Möglichkeiten des Zusatzpaketes aufgezeigt.
Wenn Sie das Zusatzpaket installiert haben, beachten Sie bitte den Bereich Initialisierung am Endes des Dokumentes.

Einführung

Wer benutzt  Wavelets ?
Jeder, der Signale nach ihrem Frequenzgehalt hin untersucht, sowie für Daten- und Bildkompression.
Hier einige Anwender:

Was sind Wavelets ?

Der Begriff Wavelet bedeutet einfach kleine Wellen oder  Wellenpakete.  Hier sind einige Beispiele :

Die Sinus-Welle setzt sich (periodisch) ins unendliche fort. Die anderen Wellen sind Wavelets und fallen schnell auf Null ab. Die Vorteile davon sind, daß man mit Hilfe von Wavelets lokale Effekte eines Signals darstellen kann, wohingegen die gewöhnliche Sinus-Welle über den gesamten Bereich ausgedehnt ist.

Herkunft und Gründe für die Entwicklung von Wavelets

Der Gebrauch von Wavelets hat eine Vielzahl von Methoden für das Analysieren, Synthetisieren, "De-Noising" und für die Kompression von Signalen und Bildern hervorgebracht oder erweitert. Wissenschaftler und Ingenieure untersuchen zeitabhängige Daten meistens in der  Frequenz-Domäne (Oszillationen). Das Verfahren dazu ist bekannt als die Fourier Transformation.  Diese jedoch ist durch  die  unendliche Ausbreitung der Sinus-Welle beschränkt. Deshalb kann die Fourier Analyse nur begrenzt mit Signalen umgehen, die abrupte Änderungen mit der Zeit haben.

Wavelets sind Wellen, die schnell auf Null abfallen und deshalb gut lokalisiert sind. Die Zerlegung einer Funktion auf diese Weise hat Gebiete wie die Bildverarbeitung, Optik, Elektronik und die Aerodynamik revolutioniert.

Hauptanwendungs-Gebiete

• Kompression

• "De-noising" (Entfernen von Rauschen)

• Zeit-Frequenz Analyse

Für diese Präsentation werden wir die mathematischen Hintergründe der Wavelets außer acht lassen. Stattdessen sollen durch einige wichtige Beispiele die Gründe für den Gebrauch von Wavelets dargestellt werden.

Fourier Analyse

Untersuchung der Frequenzen eines Signals

Mathematica hat dafür schon eine eingebaute Funktion namens Fourier.  Diese vollführt eine diskrete Fourier Transformation der Beispieldaten.

Als ein Beispiel, betrachten wir eine periodische Sinus-Welle, versehen mit etwas zufälligem Rauschen.

ListPlot[data,  PlotJoined  ->  True];

Nun erzeugen wir die Fourier-Transformierte.

Durch die Darstellung eines Leistungsspektrums kann man die vorherrschenden Frequenzen des Signals identifizieren.

ListPlot[Abs[(ft^2)],  PlotJoined  ->  True,  
      PlotRange  ->  All,  
      AxesLabel  ->  {"Frequency",  "Energy"}];

Klar erkennbar ist die Frequenz bei 30 Hz. Der zweite Peak resultiert von der periodischen Natur der Welle (256-30).

Die Fourier Analyse ist ein exzellentes Werkzeug, um den Frequenzgehalt über das gesamte Signal zu analysieren.

Darstellung eines Signals aus Cosinus und Sinus Funktionen

Jede periodische Funktion/Signal kann durch eine Überlagerung von trigonometrischen Funktionen (Sinus und Cosinus Basis) dargestellt werden. Mit Mathematica kommt ein Standard Erweiterungspacket für die Fourier Analyse:

Um die Näherung einer Funktion mittels Fourier Analyse  zu erläutern, betrachten wir  die Sprung-Funktion am Ursprung:

pa  = Plot[f[x],  {x, -4, 4}];

Schauen wir uns an, wie das Signal durch ein paar wenige Terme reproduziert wird.

fts  =  FourierTrigSeries[f[x], {x, -4, 4}, 15];

pb  =  Plot[fts,  {x, -4, 4}];

Die Kombination dieser beiden Plots zeigt das Problem.

Show[{pa, pb}];

Hier zeigt sich die Limitierung der Fourier Analyse. Neben der Tatsache, daß wir eine unendliche Zahl von trigonometrischen Termen benötigen, um das Signal exakt darzustellen, tritt auch noch ein seltsamer Effekt direkt neben dem Sprung (Diskontinuität) auf. Dieser wird als Gibbs Phänomen bezeichnet.

Beschränkungen der Fourier Analyse

• Daten die scharfe Peaks oder Diskontinuitäten besitzten können zum Gibbs Phänomen führen (siehe oben).

• Um genaue Ergebnisse zu erhalten, wird eine große Zahl von trigonometrischen Termen benötigt.

• Die Fourier Analyse hat Probleme mit lokalen Zeit-Frequenz Wechseln innerhalb des Signals.

Wavelets

Anmerkung

Möchte ein Ingenieur feststellen, welche Frequenzen einem Signal zugrunde liegen, so ist die Standard Fourier Analyse ausreichend. Ist man jedoch daran interessiert, wie sich die Frequenz mit der Zeit ändert, so sind Wavelets angebrachter.

Wavelets bieten eine lokale Frequenzabschätzung über ein endliches Zeitintervall.

Deshalb sind Wavelets gut für die Darstellung von Funktionen geeignet, die starke zeitliche Änderungen aufweisen.
Die Fourier Analyse kann solche lokalen Effekte nicht behandeln - wie schon gezeigt worden ist.

Darstellung von Signalen

Betrachten wir ein gut lokalisiertes Stoß-Signal (wie es sich beim Lösen einer Gleichung ergibt, die häufig in der Dynamik von Flüssigkeiten Verwendung findet):

spl  =  ListPlot[data,  PlotJoined  ->  True];

Zuerst müssen wir die Koeffizienten des Wavelets generieren. Dies ist einfach eine Liste. Wir können diese Liste von Koeffizienten benützen, um das tatsächliche Wavelet zu konstruieren.

filter  =  DaubechiesFilter[2]

Möchten wir einen Blick auf das hieraus gewonnene Wavelet werfen, so kann man den folgenden Befehl verwenden.

ListPlot[Evaluate[Wavelet[filter, 5]],  PlotJoined  ->  True];

Durch ändern der Nummer "5" erhält man eine stärker/weniger detaillierte Grafik.

Als nächstes benutzen wir die Koeffizienten des Wavelets, um die  Wavelet Transformation durchzuführen. Mit dem Wavelet selber brauchen wir uns dabei nicht zu beschäftigen!

Die Wavelet Transformation umfaßt eine Liste von Listen.

Length /@ wtdata

Eine Möglichkeit, die Positionen der Wavelet-Koeffizienten zu veranschaulichen, ist der Befehl PlotCoefficients.

PlotCoefficients[wtdata,  Frame  ->  True,  
      PlotRange  ->  All];

Jede Linie entspricht der Position der Wavelet-Koeffizienten, welche die Daten darstellen. Die Länge der Linien bezieht sich auf den Wert der Koeffizienten. Wie man sieht, sind die meisten Wavelets um das Zentrum des Stoßes herum angeordnet. Die unterste Zeile repräsentiert den Residuen-Trend des Signals. Die oberen Linien entsprechen den feineren Details des Signals.

Eine weitere Möglichkeit dies zu sehen stellt der Befehl PhaseSpacePlot dar. Dies ist eine Zeit-Frequenz Graphik der Wavelet Transformation.

PhaseSpacePlot[wtdata,  Frame  ->  True,  
      FrameLabel  ->  {"Time",  "Frequency"}];

Jedes Rechteck entspricht einem Koeffizienten. Je dunkler das Rechteck, desto größer der Koeffizient. Alle Rechtecke haben die gleiche Fläche. Je lokalisierter das Signal wird, desto unbestimmter wird die Frequenz. Dies ist als das Unbestimmtheits Prinzip bekannt. Beispielsweise ist die Struktur nahe des Zentrums gut lokalisiert, aber die tatsächliche Frequenz ist es nicht.

Die obige Graphik zeigt, daß ungefähr 15-20  Rechtecke wichtig sind. Wir komprimieren deshalb alle Koeffizienten außer den 20 größten:

Die Inverse Wavelet Transformation dieser 20 Koeffizienten alleine führt zu:

ListPlot[newdata,  PlotJoined  ->  True];

Das Ergebnis kommt den Originaldaten sehr nahe.

Show[{%,  spl}];

Vergleichen wir diese Kompression mit der Fourier Analyse, indem wir zuerst die diskrete Fourier Transformation anwenden.

Anwendung der gleichen Kompression

Inverse Transformation

Benutzung des Befehls Chop, um die bei der Inversen Transformation auftretenden, sehr kleinen Zahlen zu eliminieren.

ListPlot[Chop[fnew]];

Nun kann man erkennen, wie schlecht die Fourier Analyse Stoß-ähnliche Signale wiedergibt.

Show[%,  spl];

Wieder tritt bei der Fourier Analyse das Gibbs Phänomen auf.

Kompression mit Wavelets ist besser!

Entfernen von Rauschen

Als nächstes betrachten wir verrauschte Daten. Das Problem dabei ist, daß die Daten selber scharfe Peaks besitzen (ohne dem Rauschen). Hält man sich an die Standard Fourier Analyse für das Entfernen von Rauschen, so werden dabei einige der den Daten zugrundeliegenden Strukturen entfernt. Dies kann den Ingenieur zu dem Schluß bringen, daß die zugrundeliegenden Daten glatt sind.

orig  =  ListPlot[data,  PlotJoined  ->  True];

Es gibt hier einige sehr scharfe Peaks innherhalb dieser Daten. Mit dem Wavelet Explorer hat der Ingenieur die Wahl über das Glätten. Er kann entweder die zugrundeliegende Struktur beibehalten, oder sie wie die Fourier Analyse ausglätten.

Es sein nun angenommen, daß die Daten inhärent nicht-glatt sind, und wir einfach das Rauschen entfernen wollen.

Wie zuvor starten wir, indem wir einen Filter auswählen.

Anwenden der Wavelet Transformation.

Length /@ wtdata

Die erste Liste ist der Residuen-Trend der Daten. Wir möchten diejenigen Wavelets, die diesen Teil der Daten darstellen, nicht komprimieren.
Der Befehl Rest[] sagt Mathematica, die erste Liste zu ignorieren.

Die gewählte Zahl besagt, daß alle Wavelet-Koeffizienten darunter verkleinert werden. Größere Zahlen (z.B. 5) machen das Ergebnis glätter. Es sei hier daran erinnert, daß die Daten als nicht-glatt angenommen werden.

dnoi  =  ListPlot[%,  PlotJoined  ->  True];

Man beachte, wie die scharfe Struktur des Originals beibehalten wird. Die Fourier Analyse kann nur schlecht mit solch scharfen Änderungen umgehen und wird deshalb das Ergebnis überglätten.

Show[GraphicsArray[{{dnoi},  {orig}}]];

Selbstverständlich kann man die Wavelets dazu benutzen noch weitere Stufen zu glätten.

ListPlot[InverseWaveletTransform[ %,  s4],  
      PlotJoined  ->  True];

Show[GraphicsArray[{{%,  orig}}]];

Wavelet-Paket Analyse

Wavelets können dazu benutzt werden um zu klären, wie sich die Frequenz mit der Zeit ändert. Die Standard Fourier Analyse kann dies nicht, da die zeitliche (oder räumliche) Komponente des Signals im Gegensatz zu Wavelets nicht lokalisiert ist. Bis in die 90er Jahre wurde das Signal in gleichen Abständen in "Scheiben" geschnitten. Dieses Verfahren ist als die 'Windowed Fourier Transform' bekannt und unterliegt an den Grenzen der Scheiben auch den Problemen des Gibbs Phänomens.

Deshalb müssen wir uns den Wavelets zuwenden. Wir werden Wavelet-Pakete zur Kompression und für Zeit-Frequenz Experimente anwenden.

Hier analysieren wir die Kompression eines Sprachsignals welches das Wort Mathematica bildet.

ListPlay[mfile,  SampleRate  ->  8012,  PlayRange  ->  All];

Das Signal beinhaltet verschiedene "Pakete" von Daten mit unterschiedlichen Frequenzen.  Die Länge der Daten ist:

Length[mfile]

Probieren wir einen anderen Filtertyp aus. Der folgende wird gewählt, da er nahezu (aber nicht vollständig) symmetrisch ist. Ein bißchen wie das Signal, das wir betrachten.

Sobald wir Daten vorliegen haben, die Gruppen von Wellen (genannt Wellen-Pakete) beinhalten, ist die Wavelet-Paket Transformation besser als die normale Wavelet Transformation. Ohne sich im Detail zu verlieren, kann man sagen, daß die Wavelet-Paket Transformation die höheren Frequenzen innerhalb von Wellen-Paketen behandeln kann.

Hier ist die Wavelet-Paket Transformation der Daten:

Eine Kompression des Signals auf nur 10% der ursprünglichen Welle und die Anwendung der Inversen Wavelet-Paket Transformation ergibt:

ListPlay[iwp,  SampleRate  ->  8012,  PlayRange  ->  All];

Somit hat man trotz einer Komprimierung auf nur 10% der Originalgröße eine gute Darstellung des Signals beibehalten.

Eine weitere Option des Wavelet Explorers ist die Untersuchung des lokalen Frequenzgehalts eines Signals.

ListPlay[cdata,  PlayRange  ->  All,  SampleRate  ->  8012];

Dieses Signal hat verschiedene Frequenzen. Um diese so gut wie möglich zu isolieren, benützen wir die Wavelet-Paket Transformation und den PhaseSpacePlot.

Hier benutzen wir eine Cosinus Wellen-Paket und schneiden die Daten in 4 Hauptsegmente. Die Länge des Signals is 4096.

Nun nach dem Ausführen der Transformation können wir das Ergebnis in Frequenz-Zeit Form anschauen.

PhaseSpacePlot[cp,  Frame  ->  True,  ExchangeAxes  ->  True,  
      FrameLabel  ->  {"Time",  "Frequency"}];

Wiederum ist eine Unschärfe zwischen der lokalen Frequenz und der Zeit vorhanden.

Allgemeine Beispiel aus der Praxis

Seismische  Ereignisse

Das Folgende sind Daten, die ein seismisches Ereignis darstellen:

Erneut benutzen wir die Wavlet-Paket Transformation und bestimmen das Zeit-Frequenz Verhalten.

Zeit-Frequenz Ergebnisse

Hier berechnet wir das Frequenz-Zeit Verhalten unter Benutzung von Wavelets.

Show[GraphicsArray[{{ssd}, {psp}}],  
      DisplayFunction  ->  $DisplayFunction];

Die höheren Frequenzen sind, wie zu erwarten, beim Einsetzen der seismischen Erschütterung.

Entfernen des Rauschen von Daten über Sonnenflecken.

Das Folgende gibt die Sonnenfleckenaktivität der letzten 200 Jahre an.

sundata  =  ListPlot[data,  PlotJoined  ->  True];

Length[data]

Für die Wavelet Analyse müssen die Daten eine Länge von (2^n) haben. Deshalb füllen wir das Ende der Daten mit Nullen auf.

sunspots  =  Join[data, Table[0., {173}]];

Nun wählen wir einen Filter für die Wavelet Transformation aus. Wir werden mit verschiedenen Wavelet-Filtern zum Entfernen des Rauschen experimentieren. Da wir das tatsächliche Modell des Rauschen nicht kennen, ist ein experimenteller Zugang ideal.

({{DaubechiesFilter[2]; }, {DaubechiesFilter[12]; }, {LeastAsymmetricFilter[6]; }, {BiorthogonalSplineFilter[10, 4]; }})

Nun wenden wir die Wavelet Transformation an.

Die Wavelet Transformation ist eine Liste von Listen in Mathematica. Wenn wir diese auflösen, können wir die Verteilung der Wavelet-Koeffizienten betrachten.

ListPlot[Sort[Flatten[wtdata]],  
      PlotRange ->  {All,  {-150, 150}}];

Aus diesem Graphen ist eine Schwelle (zur Kompression der Koeffizienten) von ungefähr 50 als eine gute Wahl zu entnehmen.

lp  =  
      ListPlot[newdata,  PlotJoined  ->  True,  
        PlotRange  ->  {All,  {0, 250}}];

Show[GraphicsArray[{sundata, lp}]];

2-D data compression

Vor einiger Zeit stand das FBI vor dem Problem all seine gesammelten Fingerabdrücke zu digitalisieren. Dies benötigte zu viel Speicherkapazität. Zu dieser Zeit basierte das JPEG Kompressionsverfahren auf der "Windowed Fourier Transform" - Zerteilung der Daten in Schnitte und Benützung von Sinus und Cosinus Funktionen. Die hohen Kompressionsraten erreichbar durch Wavelet Transformation sind nun Standard für das FBI.

Betrachten wir einen einzigen Fingerabdruck.

ListDensityPlot[data,  Mesh  ->  False];

Wieder gibt es meherere Wahlmöglichkeiten für den Filter.

Wir werden sowohl eine Wavlet Transformation als auch eine Cosinus Transformation berechnen. Die letztere wird für den JPEG/MPEG Kompressions Algorithmus verwendet. Um den Effekt der Kompression zu zeigen, werden wir nur 7% der Originaldaten behalten!

Mit der Cosinus Transformation benützen wir Scheiben der Daten (Fenster).

Nun können wir  die Cosinus Transformation mit der Wavelet Transformation vergleichen. Wir benutzen nur 7% der Originaldaten!

ListDensityPlot[ict,  Mesh  ->  False];

Das nächste Ergebnis ist von den Wavelets.

ListDensityPlot[iwt,  Mesh  ->  False];

Wie man sehen kann, ermöglichen die Wavelets eine besseres Kompressionsverhältnis

Initialisierung

Pfadeinstellung

Dieser Pfad muß auf die Wavelet-Daten zeigen, wenn Sie das Zusatzpaket besitzen und nicht nur diese Demonstration.
Beispiel: 
"D:\\Programme\\Wolfram Research\\Mathematica\\4.0\\AddOns\\Applications\
     \\WAVELETS\\DATA"

WaveletPath  =  "D:\\Programme\\Wolfram \
Research\\Mathematica\\4.0\\AddOns\\Applications\\WAVELETS\\DATA";

Copyright

Dr. Gerard J. Cafolla
Wolfram Research Europe Ltd.
Übersetzung und Modifikation ADDITIVE GmbH
Andreas Heilemann
©2000 all rights reserved.