Hadamard Transformation

[Peter K. (peter26)]

Hallo,

Ich habe zufällig bei dem MLS - Audiomessverfahren gelesen das die 
Kreuzkorrelations mit der Hadamard Transformation realsiert werden kann 
... ich habe leider noch zu wenig Ahnung bzwl. dieser Trafo ...

Ist dies auch für normale Audiosignale (z.b zur Time Delay Estimation ) 
einsetzbar ? bzw. wenn ja gibts hier vl. gute Literatur zu dem Thema  ?

Bzw. gibts andere schnelle Kreuzkorrelationsverfahren (derzeit arbeite 
ich einfach mit der Mulitiplikation der Spektren) ?

Dankeschön
lg Peter

[Detlef _. (detlef_a)]

Just eine ähnliche Frage gabs vor kurzer zeit mal hier:
Beitrag "Schnelle Faltung"

So richtige Antworten gabs nicht. Ich bin der Meinung, dass man mit der 
Hadamard Transformation weder Korrelationen noch Faltungen rechnen kann. 
Aber sehr tief geschaut habe ich darein nicht und würde mich gerne eines 
Besseren belehren lassen. Hast Du links?

Cheers
Detlef

[peter26 (Gast)]

hi ... hab leider auch nur die hinweise das die trafo für die 
kreuzkorrelation bei der mls methode verwendet werden kann (wobei nicht 
klar ist warum ein mls signal so besonders ist ?! ... ist ja im weiteren 
sinn auch nur ein audio signal ...

lg peter

[Bernd (Gast)]

Was genau meinst du mit "Mulitiplikation der Spektren" ?
Es gibt viele Methoden die KKF zu bestimmen. Wie es scheint, hast du 
einen starken µC am laufen?

http://books.google.de/books?id=tYCJ4ZEGhsMC&lpg=PP1&dq=messtechnik&pg=PA215#v=onepage&q=&f=false

[Hagen R. (hagen)]

Ja das geht alles was du dir vorstellst. Ich selber habe das mit der 
Soundkarte schon realisiert. Zwei Lautspreicher, also Stereo werden mit 
zwei unterschiedlichen MLS gefüttert (hab das aber BPSK moduliert). Ein 
Mono Mikrophon sampelt nun mit 48KHz dieses Soundgemisch und berechnet 
live per Hadamard Transformation die Korrelation für beide MLS. Aus 
deren Daten wird die Laufzeit des Schalls für jede MLS berechnet und da 
der Abstand der beiden lautsprecher bekannt ist ortet man die Position 
des Mikrophones in der 2D Ebene. Ging auf 3mm genau.

Man kann also sehr wohl MLS per Hadamard Kreuz-/Autokorrellieren. 
Wichtig ist dabei die sogennannte Eingangs/Ausgangs-Permutation der 
Meßwerte. Man kann also nicht wie bei der FFT das kontinuierliche 
Eingangssignal, zb. 1024 Soundsamples, in die HT einspeisen sondern muß 
die Samples umsortieren. Diese Permutation, eg. Umsortierung, der 
Samples und Ausgangsdaten der HT ist das Geheimnis dabei.

Es gibt unterschiedliche Wege diese Permutationsmatirzien zu berechnen. 
Ich bevorzuge eine ziemlich einfache Idee. Zu jeder MLS existiert eine 
inverse MLS. Aus dem Polynom der MLS kannst du sehr einfach das Polynom 
dieser inversen MLS berechnen. Nun erzeugst du direkt aus dieser 
inversen MLS die Permutationsmatix.

Falls du was mit Delphi Source anfangen kannst, siehe das Attachment.

Gruß hagen

PS: ürbigens nimmt man deshalb MLS weil sie orthogonal sind und das ist 
zwingend für die FHT

[Peter K. (peter26)]

Bernd schrieb:
> Was genau meinst du mit "Mulitiplikation der Spektren" ?

Hi ... ich meine damit das ich meine Kreuzkorrelation derzeit so 
realsierit habe das ich die beiden Spektren miteinander multipliziere 
(natürlich konjugiert komplex) -> Generalized Cross Correlation.

lg

[Peter K. (peter26)]

Hagen Re schrieb:
> PS: ürbigens nimmt man deshalb MLS weil sie orthogonal sind und das ist
> zwingend für die FHT

Hi ... dankeschön für die ausführliche Antwort ...

Ich habe mich vl. ein wenig zu ungenau ausgedrückt bzgl der FHT ... ich 
möchte diese Transformation gern zur Berechnung der Kreuzkorrelation 
normaler Audiosignal verwenden (Mikrofonarray) also verwende ich keine 
MLS Signale -> ergibt sich die Kreuzkorrelation dann auch einfach aus 
dem Produkt der Spektren (also S * S' -> und dann inverse FHT)  ?

Dankeschön vielmals
Peter

[Michael O. (mischu)]

Die FHT ist meines Verständnisses nach ein Transformation im 
Zeitbereich, die durch einen Butterfly-Algorithmus eine schnelle Faltung 
des Eingangssignals realisiert, allerdings nur für binärwärtige MLS 
Folgen.
Die MLS hat ein weißes Frequenzspektrum. Alle Übertragngsfunktionen in 
der Messkette verändern dieses Spektrum. Die gemessene spaktral 
veränderte Signalfolge kann mittels der FHT in eine Impulsantwort 
umgerechnet werden und enthält den Frequenzgang aller Elemente.

Allerdings kannst Du nur MLS als Anregungssignale verwenden, damit die 
FHT die Rücktransformation der MLS in eine Impulsantwort erledigt.
Willst Du beliebige Anregungssignale nehmen, kommst Du nicht um eine FFT 
des Messignals mit anschließender Multiplikation (des Kehrwertes der 
Anregungsfunktion) herum, um die Übetragungsfunktion zu ermitteln.

http://books.google.de/books?id=yIRrs1npLRAC&pg=PA1106&lpg=PA1106&dq=fht+mls&source=bl&ots=AHYfQ9msxa&sig=Ztm9dlfloIWoz6oxefNaJYBFoE8&hl=de&ei=689vS8DsM5fsmwOr1vHKBA&sa=X&oi=book_result&ct=result&resnum=5&ved=0CBwQ6AEwBA#v=onepage&q=fht%20mls&f=true

[Detlef _. (detlef_a)]

Hallo,

ich wärme den thread nochmal auf, weil ich noch was zur Korrelation mit 
der Hadamard Transformation habe und das interessierten Kräften nicht 
vorenthalten möchte.

Man kann mit der Hadamard Transformation bedingt Korrelationen rechnen. 
Das ist für Maximum length Sequenzen (MLS) mit 
Dirac-Autokorrelationsfkt. möglich. Gezeigt wird das in zwei papern von 
Jeffrey Borish, die in dem Weinzierl-Buch zitiert werden, auf das mein 
Vorposter mischu hingewiesen hat.

Wenn man eine MLS mit Dirac-Autokorrelationsfkt. in ein System schickt 
(im vorliegenden Fall eine Raumakustik samt Lautsprecher) und dann das 
resultierende Signal mit der MLS korreliert kommt die Impulsantwort des 
Systems raus. Diese Korrelation wird mit der Hadamard Transformation 
berechnet, das geht auf schwacher Hardware.

Ich weiß nicht, ob der OP noch interessiert ist, aber ich habe 
jedenfalls viel dabeigelernt ;-)).

Angehängt das relevante Papier und meine Matlab-Basteleien dazu, die das 
Beispiel des papers nachprogrammieren.

Cheers
Detlef

[Harald M. (mare_crisium)]

Detlef,

na klar interessiert das Thema noch :-) ! Ich war in der Zwischenzeit 
auch nicht faul und habe versucht (eigentlich für den Privatbedarf, 
deshalb nicht über Fehler erhaben ;-) ) mir den Zusammenhang zwischen 
MLS- und Hadamard-Matrizen zu verdeutlichen. Das Resultat hänge ich zur 
kritischen Diskussion mal als .pdf an. - Deine Literaturstelle ist aber 
auch nicht von schlechten Eltern ;-).

Deine Hadamard-Versuche kann ich mangels MatLab leider nicht goutieren. 
Ich bin aber daran, zu versuchen, einen geschlossenen Ausdruck für die 
Koeffizienten auszubaldowern, die man bekommt, wenn man die 
Fourier-Funktionen mit Hadamard-Funktionen analysiert. Diese 
Koeffizienten sind natürlich komplex. - Wie gesagt, das ist noch sehr 
stark im Versuchsstadium (Ansatz im zweiten .pfd). Ich bin auch noch 
nicht überzeugt, dass es überhaupt geht.

Und nicht vergessen:

Math rulez!

Ciao,

mare_crisium

[Detlef _. (detlef_a)]

Hallo Harald,

>>Deine Hadamard-Versuche kann ich mangels MatLab leider nicht goutieren.

Zieh Dir Octave oder Scilab, da kommst Du für lineare Algebra nicht 
drumrum.

Dein Hadamard*.pdf geht genau in die Richtung, die Meister Borish 
einschlägt, die von Dir zitierte Diplomarbeit hat das paper auch im 
Literaturverzeichnis, es werden diegleichen Beispiele verwendet. Das 
paper schien mir die Dinge ganz gut zu erklären, zumindest hab ichs 
soweit verstanden, dass ich es programmieren konnte.

Da ich keine Raumakustiken vermessen will ist mein Erkenntnisinteresse 
erstmal gestillt.

>>Und nicht vergessen:
>>Math rulez!

yep :-))))

Cheers
Detlef

[Detlef _. (detlef_a)]

Hallo,

ich möchte den thread nochmal aufwärmen, weil ich mit der Hardamard 
Transformation und MLS Sequenzen mal ernsthaft Impulsantworten berechnen 
will.

Mit dem Matlab script, das ich in diesem thread vor mehr als vier Jahren 
:-/ gepostet habe geht das auch im Prinzip, es sind jedoch sehr große 
Matrizen notwendig, für eine 13-Bit MLS der Größe 8191x8191, für eine 
Produktivversion auf schwacher hardware is das nix.

Diese Matrizen sortieren den Meßvektor nur um, so wie Hagen Re das auch 
geschildert hat.

Also: Hagen, kannst Du noch mal genauer darlegen, wie Du die 
Umsortierung mit Deiner 'inversen MLS' machst. Was meinst Du mit 
'inverser MLS', hast Du Literaturzeiger.

Ich hoffe, Du bist hier noch aktiv.

THX
Cheers
Detlef

[Detlef _. (detlef_a)]

Ich habe die Hadamard Transformation nochmal angewandt um mit einer 
Soundkarte ein Mikrophon zu orten. Hatte Hagen schon so gemacht, klappt 
gut bis in den cm-Bereich.

Die beiden Lautsprecher senden jeweils eine MLS. Das Mic nimmt beide 
Signale auf und ich bestimme mit Hadamard die Impulsantwort der beiden 
Kanäle. Im angehängten Bild sieht man die beiden Impulsantworten. Aus 
den Maxima der Impulsantworten und der Geometrie der Anordnung läßt sich 
die Position des Mic berechnen.

Ich verwende MLS der Länge 511 Abtastwerte. Vor dem Senden filtere ich 
die noch mit einem Hochpaß von ca. 5kHz Grenzfrequenz damit's in der 
Wohnung nicht so laut wird ;-) . Das Empfangssignal muß für beide MLS 
vorsortiert werden, wie Hagen das schon gesagt hatte. Dann geht es durch 
eine FHT (fast hadamard transformation) und dann muß noch mit einer 
Matrix nachmultipliziert werden. Heraus kommen letztlich die 
Impulsantworten der beiden Kanäle. Die Herleitung der Vorsortierung und 
der Nachmultiplikation ist extrem speicheraufwendig, bei 511er Länge war 
bei meinem Rechner Schluß. Wg. SNR und Robustheit der Ortung sind aber 
längere MLS wünschenswert. Das SNR ist nicht so berauschend, das Mic hat 
auch eine zweifelhafte Qualität und eine ungünstige Richteigenschaft.

Ich schreibe das hier, weil ich das für ein extrem interessantes 
Verfahren zur Positionsbestimmung halte. Es ist in der Durchführung 
nicht aufwendig und liefert gute stabile Ergebnisse.
Damit kann man interessantes Zeug basteln, gibt vllt. auch nen gutes 
Thema für ne Bachelor-/Master Arbeit ab.

War Spaß gewesen.

Cheers
Detlef

[Harald M. (mare_crisium)]

Detlef _a schrieb:
> ... und dann muß noch mit einer
> Matrix nachmultipliziert werden...

Hallo, Detlef,
freut mich, dass das Verfahren bei Dir jetzt klappt :-). Aber was für 
eine Matrix meinst Du? Die Permutationsmatrix zum Umordnen des 
Ergebnisvektors? Das ist doch auch nur eine Permutation und lässt sich 
sequentiell (nicht als Matrixmultiplikation) bequem ausführen.

Sendest Du auf beiden Kanälen diesselbe MLS?

Ich verwende das Verfahren bisher mit nur einem Kanal als einfachen 
Entfernungsmesser/ Hinderniswarner für ein Asuro-artiges Fahrzeug. 
Sender und Empfänger sind bei mir 40kHz-Ultraschallkapseln. Weil mich 
nur die Amplitude interessiert (nicht die Phasenlage), werte ich das 
verstärkte, vollweggleichgerichtete Empfangssignal aus. Das Ganze läuft 
in Assembler mit einer 255-Bit MLS (wg. Speicherbedarfs) klaglos auf 
einem ATmega644 bei 16MHz.

Ciao, Harald

P.S.: Habe mich ausserdem nochmal mit der Frage beschäftigt, ob der 
Faltungssatz der Fouruertransformation auch für die 
Hadamardtransformation gilt. Scheint nicht so zu sein. Habe dabei 
allerdings den reverse bit shuffle nicht berücksichtigt. Ich hatte 
eigentlich gemeint, der Faltungssatz gelte für alle vollständigen, 
orthogonalen Funktionssysteme.
Ha.

[Detlef _. (detlef_a)]

Hallo Harald,

Harald M. schrieb:
> Detlef _a schrieb:
>> ... und dann muß noch mit einer
>> Matrix nachmultipliziert werden...
>
> Hallo, Detlef,
> freut mich, dass das Verfahren bei Dir jetzt klappt :-). Aber was für
> eine Matrix meinst Du? Die Permutationsmatrix zum Umordnen des
> Ergebnisvektors? Das ist doch auch nur eine Permutation und lässt sich
> sequentiell (nicht als Matrixmultiplikation) bequem ausführen.

Ja, der Ergebnisvektor muss bei mir nochmal mit einer Matrix 
multipliziert werden. Ich habe probiert, das auch auf eine Permutation 
der Vektorelemente zurückzuführen, das ist mir allerdings trotz 
allerschärfster Bemühungen nicht gelungen.

Details, entsprechend obigem Matlab script:
Ich berechne die Impulsantwort hr aus der Messung y mit folgender 
Formel:
hr=inv(MDnn)*P2*S2*H*S1*P1*y
Das
S2*H*S1*P1*y
ist dabei kein Problem: P1 ist eine Permutation der Elemente von y, S1 
wie auch S2 ist nichts Wesentliches, H ist die Fast Hadamard 
Transformation.

inv(MDnn)*P2 hat es aber in sich. P2 berechnet sich so:
bla=N *P1.';
P2=(G*inv(bla)).';

N ist die Matrix der rotierten MLS Folgen, P1 ist die schon bekannte 
Permutation der Eingangselemente. Aber was ist dann P1.', die inverse 
Permutation wäre inv(P1) ???
G ist die Hadmard Transformation ohne die erste Spalte und die erste 
Zeile, aber was ist P2=(G*inv(bla)).'? Und dann kommt noch inv(MDnn), 
das ist die inverse der Matrix der rotierten Autokorrelationsfkt.., wie 
krieg ich die denn aufgelöst?

P2 sieht auch komisch aus, hat paar +/-1 und positive und negative 
Kehrwerte von Zweierpotenzen. Riecht nach Vereinfachbarkeit, ist meinem 
Intellekt allerdings nicht gelungen :-/ .
>
> Sendest Du auf beiden Kanälen diesselbe MLS?
>

Nein, zwei verschiedene. Dann kann ich für beide Kanäle die 
Impulsantwort getrennt berechnen und die Entfernung zu BEIDEN 
Lautsprechern berechnen, i.e. eine Ortung in der Ebene machen.

> Ich verwende das Verfahren bisher mit nur einem Kanal als einfachen
> Entfernungsmesser/ Hinderniswarner für ein Asuro-artiges Fahrzeug.
> Sender und Empfänger sind bei mir 40kHz-Ultraschallkapseln. Weil mich

Wie geht das, die Kapseln sind total schmalbandig, die sind bei 39kHz 
schon 20dB im Keller !? Die MLS ist breitbandig. Modulierst Du einen 
40kHz Träger langsam mit der MLS?

> nur die Amplitude interessiert (nicht die Phasenlage), werte ich das
> verstärkte, vollweggleichgerichtete Empfangssignal aus. Das Ganze läuft
> in Assembler mit einer 255-Bit MLS (wg. Speicherbedarfs) klaglos auf
> einem ATmega644 bei 16MHz.
>
> Ciao, Harald

Verstehe ich komplett nicht?? Je länger die MLS ist, umso besser wird 
das SNR. BTW.: Ich habe hier noch breitbandige Ultraschall Speaker/Mics 
rumliegen, die ich noch nicht ausprobieren konnte 
http://www.knowles.com/eng/Products/Microphones/Surface-mount-MEMs
Wunderbare Technik: klein, rein digital (die machen als Ausgangssignal 
eine pulse density modulierte Folge, wie ein sigma/delta Wandler), 
breitbandig 100kHz, billig (1E oder so). Mit denen kann man schöne 
Fledermausbelauscher und MLS Entfernungsmesser bauen.

>
> P.S.: Habe mich ausserdem nochmal mit der Frage beschäftigt, ob der
> Faltungssatz der Fouruertransformation auch für die
> Hadamardtransformation gilt. Scheint nicht so zu sein. Habe dabei
> allerdings den reverse bit shuffle nicht berücksichtigt. Ich hatte
> eigentlich gemeint, der Faltungssatz gelte für alle vollständigen,
> orthogonalen Funktionssysteme.
> Ha.

Ja, ich habe in Deinen Ansatz mal reingekuckt. Die Welt der 
Transformationen ist weit und groß. Leider bin ich da mit meiner 
Mathekompetenz relativ schnell am Ende, vllt. sollte ich das wirklich 
nochmal ändern, denn:

math rulez! :-))
Cheers
Detlef

[Harald M. (mare_crisium)]

Detlef,
ich hänge hier ein älteres .pdf an, in dem ich beschrieben habe, wie ich 
das mit den schmalbandigen Ultraschallwandlern mache. Im Grunde rutsche 
ich mit der MLS-Impulsfolge über die Frequenzkennlinie der Wandler hin 
und her, dadurch ergibt sich automatisch auch eine Amplitudenmodulation. 
Nicht sehr effektiv, aber es funktioniert. - Den Vorverstärker vor der 
aktiven Vollweggleichrichtung müsste ich nochmal überarbeiten, um mit 
grösserer Verstärkung eine grössere Reichweite hinzukriegen; so wie's 
jetzt ist, schneidet der Gleichrichter prozentual zuviel vom Signal weg.

Die Umformung des hadamard-Ergebnisvektors machst Du viel zu 
kompliziert. Das hängt mit der Art zusammen, wie die Permutationen 
festgelegt werden.

Das Verfahren beruht auf der Idee, die Multiplikation des Messvektors 
mit der Matrix der rotierten MLS (MLS-Matrix) durch die schnelle 
Hadamardtransformation zu ersetzen. Das geht, weil man sowohl die 
Hadamard- als auch die MLS-Matrix in ein Produkt aus zwei rechteckigen 
Matrizen zerlegen kann. Die Zeilen der ersten und die Spalten der 
zweiten Faktormatrizen kann man durch Vertauschen in einander 
überführen.
Das ist der Grund, warum man den Messvektor und den Ergebnisvektor 
permutieren muss. Und zwar den Messvektor so, dass er in die für die 
HadamardTransformation erforderliche Form überführt wird, den 
Ergebnisvektor so, dass er von der Hadamard-Form in die 
Messergebnis-Form zurückkommt. Daher kommt wahrscheinlich das "inv" in 
Deiner Rückpermutation.

Ich hoffe, das ist wenigstens teilweise verständlich :-/

Ciao, laissez math ruler! Harald

[Tobias P. (hubertus)]

Hallo Harald und Detlef,

ich erlaube mir mal, mich hier einzuklinken :-)
Mich interessiert das Thema auch sehr. Weiter oben wurde geschrieben, 
dass man mit der MLS und Kreuzkorrelation die Impulsantwort von einem 
System bestimmen kann. Damit befasse ich mich zurzeit auch sehr intensiv 
(also mit dem Finden der Impulsantwort). Habe einiges schon gemacht, 
unter anderem Fouriertransformation auch angewandt, aber die Ergebnisse 
haben nicht sehr überzeugt. Daher wäre ich am Verfahren mit der 
Korrelation interessiert.
Mir stellt sich jetzt die Frage, wie man die MLS auf den Lautsprecher 
ausgibt? @Detlef, wie hast du das gemacht.
Damit wir uns recht verstehen: mit der MLS meint ihr eine PRBS Sequenz 
maximaler Länge, ja? Maximum Length Sequence ist mir da ein Begriff. 
Wenn ihr was anderes meint wäre ich an einer Erklärung interessiert ;-)

Eure PDFs führe ich mir jetzt noch zu Gemüte. Da steht sicher auch was 
interessantes drin.


>> math rulez! :-))
> laissez math ruler!

dito! :D

[Detlef _. (detlef_a)]

Hallo,

jetzt gehts' bei mir auch, ich bekomme für die Berechnung einfache 
Permutationen, keine Matrizen mehr. Bei dem Borish-paper von oben gibt 
es Klippen, der Algo funktioniert unten gewissen Startbedingungen der 
MLS nicht, das hat mich zwei Tage Nachdenken gekostet :-)))).

Die wesentlichen Dinge habe ich zur geneigten Ansicht als Matlab-File 
angehängt. Das script ist meine Baustelle, das ist weder schön noch 
sauber, gerne aber Nachfrage wenn Interesse besteht.

Mein Speicherausbau schafft es für dieses script bis zu MLS der Länge 
511, mein Ziel sind Ortungen mit MLS der Länge 2^17-1, dazu müßte ich ~ 
8*2^34 Byte Speicher für eine Matrix haben, das klappt so nicht ;-))), 
deswegen werde ich das entweder in Matlab oder C so umbauen, dass ich 
die Permutationen für 2^17-1 berechnen kann.

>>>ich erlaube mir mal, mich hier einzuklinken :-)
Das ist öffentlich, hier kann jeder was sagen ;-)

>>>>Habe einiges schon gemacht,
unter anderem Fouriertransformation auch angewandt, aber die Ergebnisse
haben nicht sehr überzeugt.

been there, done that, geht super. Theoretisch einfacher als Hadamard, 
dafür rechenaufwendiger, man muss multiplizieren, das geht bei Hadamard 
ohne.

>>>Damit wir uns recht verstehen: mit der MLS meint ihr eine PRBS Sequenz
maximaler Länge, ja? Maximum Length Sequence ist mir da ein Begriff.
Wenn ihr was anderes meint wäre ich an einer Erklärung interessiert ;-)

MLS steht für Maximum Length Sequence, das ist eine PRBS (pseudo random 
binary sequence?) maximaler Länge, für n taps 2^n-1. Die Funktion prnum 
in angehängtem script liefert die Dir für gegebene Bitlänge.

>>>>Mir stellt sich jetzt die Frage, wie man die MLS auf den Lautsprecher
ausgibt? @Detlef, wie hast du das gemacht.

Die PRBS sind Nullen und Einsen. Die lügst Du zu -1/+1 um, das macht der 
Algo von Borish mehrmals, so werden die in ein physikalisches 
Rechtecksignal mutiert, das gebe ich einfach aus, hört sich an wie 
Rauschen, synchron nehme ich das Mic Signal auf. Ich mach das in Matlab, 
Du kannst aber auch ein .wav File erzeugen und mit beliebigem Programm 
ausgeben/einlesen. Für die 131071 samples PRBS peile ich 96kHz 
samplefrequnz an, sodass die in einer gutes Sekunde durch sind.

>>Daher kommt wahrscheinlich das "inv" in  Deiner Rückpermutation.

Das "inv" habe verstanden, Borish erkärt das auch: Die MLS haben eine 
Autokorrelation mit peak bei 2^n-1, der Rest ist -1 . Das ist kein 
Dirac, und das inv kann das kompensieren. Für DC-freie Signale großer 
Länge spielt das aber keine Rolle, das kann man weglassen, untersuche 
ich in angehängtem script auch, wie sich das unterscheidet.

Yo, Hadamard geht, schönes Erfolgserlebnis.

math rulez! :-))))
Cheers
Detlef

[Tobias P. (hubertus)]

Hallo Detlef,

> been there, done that, geht super. Theoretisch einfacher als Hadamard,
> dafür rechenaufwendiger, man muss multiplizieren, das geht bei Hadamard
> ohne.

Ja, also von der Nachvollziehbarkeit, was da abgeht, ist Fourier 
definitiv einfacher ;-) Ich habe vom Eingangs- und Ausgangssignal die 
FFT berechnet, dann die Division berechnet und die inverse FFT 
gerechnet. Dann habe ich sowas wie ein Bode-Diagramm erhalten. 
Interessanterweise habe ich hier: 
https://www.spsc.tugraz.at/sites/default/files/PA_Mikula_mehrkanalige_Impulsantworten.pdf 
ein PDF gefunden, wo das auch gemacht wird, dessen Autor meint aber, 
dass man da vor der ifft() noch was spiegeln sollte. Was mir nicht klar 
ist wieso. Weisst du weshalb?

Der Rechenaufwand um an meine Impulsantwort oder Bode zu kommen, ist mir 
ziemlich egal, es soll nur reproduzierbar und einigermassen genau sein.

Mit der MLS sollte das ja auch funktionieren und ich würde prinzipiell 
sogar noch ein wenig bessere Ergebnisse erwarten als mit dem 
Sinus-Sweep, welchen ich bis jetzt benutzt habe.

> MLS steht für Maximum Length Sequence, das ist eine PRBS

gut, dann haben wir vom Selben gesprochen ;-)





> Die PRBS sind Nullen und Einsen. Die lügst Du zu -1/+1 um, das macht der
> Algo von Borish mehrmals, so werden die in ein physikalisches
> Rechtecksignal mutiert, das gebe ich einfach aus, hört sich an wie
> Rauschen, synchron nehme ich das Mic Signal auf.

Hast du damit eine zufrieden stellende Impulsantwort berechnen können?
Ich sehe in den Posts oben gerade nicht genau, wie du das gemacht hast. 
Oder stünde das in dem eingescannten PDF weiter oben?

Danke & Gruss Tobias.



> math rulez! :-))))
insbesondere auch mit Matlab & Co. ;-)

[Hagen R. (hagen)]

Detlef _a schrieb:
> Also: Hagen, kannst Du noch mal genauer darlegen, wie Du die
> Umsortierung mit Deiner 'inversen MLS' machst. Was meinst Du mit
> 'inverser MLS', hast Du Literaturzeiger.

Keine Literatur und ich habe wahrlich danach gesucht. Also entweder habe 
ich da was gefunden was noch nicht beschrieben wurde oder aber ich finde 
den Link in die korrejte Mathematik nicht.

Zur Erzeugung deiner MLS benutzt du ja ein Generatorpolynom. Dieses 
kannst du ja als Binärzahl ausdrücken, eben exakt so wie es bei dem 
Polynom einer CRC gemacht wird. Drehe nun die Reihenfolge der Bits um 
und erzeugt damit das "inverse Polynom" so wie auch bei CRCs bekannt. Es 
entsteht ein neues Polynom das die gleichen LFSR-Eigenschaften wie das 
Originale hat. Den Output dieses LFSRs kannst du 1 zu 1 benutzten um die 
nötigen Tabellen zur Indizierung für die Hadamard Transformation zu 
erzeugen. Da der Tabellenzugriff i.A. aber viel langsammer ist als die 
Live-Erzeugung dieser Sequenz benutzte ich diese also direkt beim 
Umsortierung vor und nach der Hadamard Tranfsormation. Hier mal mein 
Delphi-PASCAL Source.

Input ist ein Array of Integer, also ein dynamisches Array aus 
vorzeichenbehafteten Ganzzahlen.
Output das gleiche Spiel als Returnwert.
Poly das Generatorpolynom als Binärzahl deiner MLS und Seed der 
Startwert des LFSRs.

Ich habe über diesen "Trick" wirklich nichts im WEB gefunden, entweder 
weil er unbekannt ist oder aber ich nicht die richtigen mathem. 
Suchbegriffe kenne. Auf alle Fälle ist dieser "Trick" wesentlich 
eleganter als all die Methoden die ich bisher finden konnte bei denen 
man aufwendige Matrizenrechnungen durchführen muß oder sogar 
Suchalgorithmen anwendet.
Vorteil ist eben das man "On the Fly" die Indizierungsmatrizen berechen 
kann und somit ohne Vorausberechnungen auskommt und eben weniger 
Speicher benötigt.

Ich hoffe du kommt mit PASCAL klar.

function DoFHT(const Data: TIntegerArray; const Poly, Seed: Cardinal): TIntegerArray;

  function Log2(const Value: Cardinal): Integer;
  // Result := Ln2(Value)
  asm
         BSR   EAX,EAX
         JNZ   @@1
         DEC   EAX
  @@1:
  end;

  function Permute(const Data: TIntegerArray; Poly, Seed: Cardinal): TIntegerArray;
  // nutze inverse MLS zum Polynom "Poly" um die Eingangssamples umzusortieren für die nachfolgende FHT
  // "Poly" ist also das Polynom der MLS kodierten Daten
  var
    I,Bits,Size,DC,Value: Integer;
    Mask,Idx: Cardinal;
  begin
    Bits := Log2(Poly);     //   7
    Mask := 1 shl Bits;
    Size := Mask -1;        //   127
    Seed := Seed and Size;  //   $7F, 127
    Idx  := 0;
    Poly := Poly shr 1;
    Mask := Mask shr 1;     //   64
    for I := 0 to Bits -1 do
    begin
      Idx := Idx shr 1;
      if Odd(Seed) then
      begin
        Idx  := Idx or Mask;
        Seed := (Seed shr 1) xor Poly;
      end else Seed := Seed shr 1;
    end;
    SetLength(Result, Size +1);
    DC := 0;
    for I := 0 to Size -1 do
    begin
      Idx := Idx shr 1;
      if Odd(Seed) then
      begin
        Idx  := Idx or Mask;
        Seed := (Seed shr 1) xor Poly;
      end else Seed := Seed shr 1;
      Value := Data[I];
      Dec(DC, Value);
      Result[Idx] := Value;
//      Write(Idx, ',');
    end;
 //   WriteLn;
    Result[0] := 0;//DC;
  end;

  function InvPermute(const Data: TIntegerArray; Poly, Seed: Cardinal): TIntegerArray;
  // inverse Permutation der Ausgangsdaten der FHT
  var
    I,Bits,Size: Integer;
  begin
    Bits := Log2(Poly);
    Size := 1 shl Bits -1;
    Seed := Poly;
    Poly := 0;
// inverse Polynom berechenen, Seed=Original Polynom, Poly=invertirtes Polynom
    for I := 0 to Bits -1 do
    begin
      Inc(Poly, Poly + Seed and 1);
      Seed := Seed shr 1;
    end;
    SetLength(Result, Size);
    for I := 0 to Size -1 do
    begin
      if Odd(Seed) then Seed := (Seed shr 1) xor Poly
        else Seed := Seed shr 1;
      Result[I] := Data[Seed];
//      Write(Seed, ',');
    end;
//    WriteLn;
  end;

  function Butterfly(Data: TIntegerArray): TIntegerArray;
  // FHT selber
  var
    I,J,K1,K2,T,S,X,Y: Integer;
  begin
    K1 := 1;
    K2 := Length(Data) div 2;
    repeat
      Y := 0;
      OpCount := OpCount + 1;
      for I := 0 to K1 -1 do
      begin
        X := Y;
        Y := Y + K2;
        OpCount := OpCount + 4;
        for J := 0 to K2 -1 do
        begin
          Write(X:5, ',', Y:5, '; ');
          S       := Data[X];
          T       := Data[Y];
          Data[X] := S + T;
          Data[Y] := S - T;
          Inc(X);
          Inc(Y);
          OpCount := OpCount + 20;
        end;
        WriteLn;
      end;
      WriteLn;
      K1 := K1 * 2;
      K2 := K2 div 2;
    until K2 = 0;
    Result := Data;
  end;

begin
  OpCount := 0;
  Result := InvPermute(Butterfly(Permute(Data, Poly, Seed)), Poly, Seed);
end;

[Hagen R. (hagen)]

hier noch par Polynome die du benutzen kannst. P12 zB. enthält alle 
nicht reduzierbaren 12 Bit Polynome, also GF(2^12).

const
  P3: array[0..1]  of Word =
    ($0B, $0D);

  P4: array[0..1]  of Word =
    ($13, $19);

  P5: array[0..5]  of Word =
    ($25, $29, $2F, $37, $3B, $3D);

  P6: array[0..5]  of Word =
    ($043, $05B, $061, $067, $06D, $073);

  P7: array[0..17] of Word =
    ($083, $089, $08F, $091, $09D, $0A7, $0AB, $0B9,
     $0BF, $0C1, $0CB, $0D3, $0D5, $0E5, $0EF, $0F1,
     $0F7, $0FD);

  P8: array[0..15] of Word =
    ($11D, $12B, $12D, $14D, $15F, $163, $165, $169,
     $171, $187, $18D, $1A9, $1C3, $1CF, $1E7, $1F5);

  P9: array[0..47] of Word =
    ($211, $21B, $221, $22D, $233, $259, $25F, $269,
     $26F, $277, $27D, $287, $295, $2A3, $2A5, $2AF,
     $2B7, $2BD, $2CF, $2D1, $2DB, $2F5, $2F9, $313,
     $315, $31F, $323, $331, $33B, $34F, $35B, $361,
     $36B, $36D, $373, $37F, $385, $38F, $3B5, $3B9,
     $3C7, $3CB, $3CD, $3D5, $3D9, $3E3, $3E9, $3FB);

  P10: array[0..59] of Word =
    ($409, $41B, $427, $42D, $465, $46F, $481, $48B,
     $4C5, $4D7, $4E7, $4F3, $4FF, $50D, $519, $523,
     $531, $53D, $543, $557, $56B, $585, $58F, $597,
     $5A1, $5C7, $5E5, $5F7, $5FB, $613, $615, $625,
     $637, $643, $64F, $65B, $679, $67F, $689, $6B5,
     $6C1, $6D3, $6DF, $6FD, $717, $71D, $721, $739,
     $747, $74D, $755, $759, $763, $77D, $78D, $793,
     $7B1, $7DB, $7F3, $7F9);

   P11: array[0..175] of Word =
    ($805, $817, $82B, $82D, $847, $863, $865, $871,
     $87B, $88D, $895, $89F, $8A9, $8B1, $8CF, $8D1,
     $8E1, $8E7, $8EB, $8F5, $90D, $913, $925, $929,
     $93B, $93D, $945, $949, $951, $95B, $973, $975,
     $97F, $983, $98F, $9AB, $9AD, $9B9, $9C7, $9D9,
     $9E5, $9F7, $A01, $A07, $A13, $A15, $A29, $A49,
     $A61, $A6D, $A79, $A7F, $A85, $A91, $A9D, $AA7,
     $AAB, $AB3, $AB5, $AD5, $ADF, $AE9, $AEF, $AF1,
     $AFB, $B03, $B09, $B11, $B33, $B3F, $B41, $B4B,
     $B59, $B5F, $B65, $B6F, $B7D, $B87, $B8B, $B93,
     $B95, $BAF, $BB7, $BBD, $BC9, $BDB, $BDD, $BE7,
     $BED, $C0B, $C0D, $C19, $C1F, $C57, $C61, $C6B,
     $C73, $C85, $C89, $C97, $C9B, $C9D, $CB3, $CBF,
     $CC7, $CCD, $CD3, $CD5, $CE3, $CE9, $CF7, $D03,
     $D0F, $D1D, $D27, $D2D, $D41, $D47, $D55, $D59,
     $D63, $D6F, $D71, $D93, $D9F, $DA9, $DBB, $DBD,
     $DC9, $DD7, $DDB, $DE1, $DE7, $DF5, $E05, $E1D,
     $E21, $E27, $E2B, $E33, $E39, $E47, $E4B, $E55,
     $E5F, $E71, $E7B, $E7D, $E81, $E93, $E9F, $EA3,
     $EBB, $ECF, $EDD, $EF3, $EF9, $F0B, $F19, $F31,
     $F37, $F5D, $F6B, $F6D, $F75, $F83, $F91, $F97,
     $F9B, $FA7, $FAD, $FB5, $FCD, $FD3, $FE5, $FE9);

   P12: array[0..143] of Word =
    ($1053, $1069, $107B, $107D, $1099, $10D1, $10EB, $1107,
     $111F, $1123, $113B, $114F, $1157, $1161, $116B, $1185,
     $11B3, $11D9, $11DF, $120D, $1237, $123D, $1267, $1273,
     $127F, $12B9, $12C1, $12CB, $130F, $131D, $1321, $1339,
     $133F, $134D, $1371, $1399, $13A3, $13A9, $1407, $1431,
     $1437, $144F, $145D, $1467, $1475, $14A7, $14AD, $14D3,
     $150F, $151D, $154D, $1593, $15C5, $15D7, $15DD, $15EB,
     $1609, $1647, $1655, $1659, $16A5, $16BD, $1715, $1719,
     $1743, $1745, $1775, $1789, $17AD, $17B3, $17BF, $17C1,
     $1857, $185D, $1891, $1897, $18B9, $18EF, $191B, $1935,
     $1941, $1965, $197B, $198B, $19B1, $19BD, $19C9, $19CF,
     $19E7, $1A1B, $1A2B, $1A33, $1A69, $1A8B, $1AD1, $1AE1,
     $1AF5, $1B0B, $1B13, $1B1F, $1B57, $1B91, $1BA7, $1BBF,
     $1BC1, $1BD3, $1C05, $1C11, $1C17, $1C27, $1C4D, $1C87,
     $1C9F, $1CA5, $1CBB, $1CC5, $1CC9, $1CCF, $1CF3, $1D07,
     $1D23, $1D43, $1D51, $1D5B, $1D75, $1D85, $1D89, $1E15,
     $1E19, $1E2F, $1E45, $1E51, $1E67, $1E73, $1E8F, $1EE3,
     $1F11, $1F1B, $1F27, $1F71, $1F99, $1FBB, $1FBD, $1FC9);

[Detlef _. (detlef_a)]

Hallo alle Interessierten,

Hagen, schön dass Du noch aktiv bist. Deinen Code kenne ich, den hattest 
Du in 2010 ja schon gepostet und ich habe aus Deiner Delphi 
Hadamard-Butterfly Implementation meine Matlab Version gebastelt. Die 
'Permute' und 'InvPermute' habe ich nach C übersetzt aber noch nicht 
getestet oder verstanden. Inc() und Cardinal kannte ich nicht, aber is 
ja nicht so schwer.

>Auf alle Fälle ist dieser "Trick" wesentlich
>eleganter als all die Methoden die ich bisher finden konnte bei denen
>man aufwendige Matrizenrechnungen durchführen muß oder sogar
>Suchalgorithmen anwendet.

Ja, auf jeden Fall. Super Verfahren, das kannte Borish, der Author des 
von mir benutzten papers nicht. Der hat vorberechnet, so wie ich das 
auch noch mache.

>>und eben weniger Speicher benötigt.
Ja, ich benutze Folgen der Länge 2^17-1, die beiden Umsortiertabellen 
nehmen 1MByte Speicher.

Alle Polynome hatte ich für Folgen der Länge 2^13-1 auch berechnet, 
davon gibt es 630. Die nächste Woche ist mein ZweitPC noch mit der 
Berechnung aller Polynome für 2^17-1 beschäftigt, zwei kenne ich schon, 
die benutze ich für die beiden Lautsprecher.

Dein Verfahren werde ich mir kommende Woche ansehen, auch im Hinblick 
auf die Verbesserung meiner Rechenkünste bzgl. Galois-Feldern.

@ Tobias
>Hast du damit eine zufrieden stellende Impulsantwort berechnen können?
>Ich sehe in den Posts oben gerade nicht genau, wie du das gemacht hast.
>Oder stünde das in dem eingescannten PDF weiter oben?

Ja, sehr zufriedenstellend, das klappt 1a. Ich habe alles 
zusammengeschrieben und als PDF mit Sourcen mal angehängt. In dem 
eingescannten PDF des papers von Borish steht der Algorithmus so wie ich 
ihn implementiert habe.

>wo das auch gemacht wird, dessen Autor meint aber,
>dass man da vor der ifft() noch was spiegeln sollte. Was mir nicht klar
>ist wieso. Weisst du weshalb?

Hab 'spiegeln' in der Arbeit nicht gefunden, welche Seite meinst Du? 
Vllt. ist die konjugiert complexe Ergänzung gemeint, die Du vor der IFFT 
machen must um ein rein reeles Zeitsignal zu kriegen?

>Der Rechenaufwand um an meine Impulsantwort oder Bode zu kommen, ist mir
>ziemlich egal, es soll nur reproduzierbar und einigermassen genau sein.

Ist genau und reproduzierbar. Selbst mein Versuch mit dem BilligMic und 
den verrauschten Billigboxen bringt sehr gut reproduzierbare Ergebnisse.

>Mit der MLS sollte das ja auch funktionieren und ich würde prinzipiell
>sogar noch ein wenig bessere Ergebnisse erwarten als mit dem
>Sinus-Sweep, welchen ich bis jetzt benutzt habe.

Ja, definitiv. Sinus-Sweep taugt nicht viel. Ich helfe Dir gerne aufs 
Pferd dabei, das ist im Blick zurück nicht so schwierig. Die Sourcen 
sind ja in meinem hier angehängten PDF dokumentiert, kann da gerne auch 
noch expliziter werden.

War großer Spass.
math rulez!
Gute Nacht
Detlef

[Tobias P. (hubertus)]

Hallo Detlef,

danke für das PDF, muss ich mir heute Abend zu Gemüte führen.

Wegen der inversen MLS. Dazu kann ich folgendes sagen, so ein bisschen 
kommt mir das bekannt vor. Und zwar kann man ja die MLS in Form eines 
Polynoms beschreiben. Ähnlich wie bei der CRC, wo es ja auch ein CRC 
Polynom gibt.

Ich habs jetzt nicht ausprobiert, aber ich vermute, wenn man dieses 
Polynom invertiert, dann bekommt man daraus eben das Polynom der 
inversen MLS. Bei der z-Transformation braucht man manchmal auch solche 
inversen Polynome. Beispiel hierfür:

Das Polynom sei

$$1 + 0.78z^{-1} + 0.3z^{-2}$$

dann heisst das inverse Polynom dazu

$$z^{-2} + 0.78z^{-1} + 0.3$$

das ist zwar aus einem Beispiel von einer z-Trafo entnommen, aber du 
kannst das ja leicht ummünzen auf positive Exponenten, denke ich.
Es ist also nicht eine Kehrwertbildung, wie man vermuten würde (weil es 
"invers" heisst ;-) ).

> Sinus-Sweep taugt nicht viel.
naja, ich habe gestern damit noch rum gepröbelt, und einigermassen 
brauchbare Resultate erzielt, aber mit der Fourier-Methode. Problem 
dabei ist allerdings, dass das Messrauschen einen Einfluss auf die 
Resultate hat, was bei anderen Verfahren nicht so sehr ins Gewicht 
fällt.


Gruss Tobias

[Hagen R. (hagen)]

Detlef _a schrieb:
> Alle Polynome hatte ich für Folgen der Länge 2^13-1 auch berechnet,
> davon gibt es 630. Die nächste Woche ist mein ZweitPC noch mit der
> Berechnung aller Polynome für 2^17-1 beschäftigt, zwei kenne ich schon,
> die benutze ich für die beiden Lautsprecher.

Wochen beschäftigt ? Im Attachment alle 7710 primitiven Polynome für 
Degree 17. Die Berechnung hat bei mir 1,3 Sekunden gedauert ;) Alle 
276480 Polynome mit zB. Degree 24 zu finden dauert dann schon 9,7 
Minuten.

Wenn du meine Polynome verwendest dann achte darauf das sie in binär die 
vollständige Schreibweise/Codierung benutzen. Im Allgemeinen musst du 
also vor der Benutzung mit den üblichen LFSR/CRC Verfahren sie um 1 Bit 
rechtsschieben eg. durch 2 ohne Rundung teilen.

Falls du mit Delphi zurecht kommst könnte ich dir die Sourcen zur 
Berechnung solcher Polynome mailen. Die Sourcen arbeiten bis 2^31-1 also 
alles Integer Arithmetik.

Cardinal = unsigned 32bit Integer

Gruß hagen

[Hagen R. (hagen)]

übrigens schönes PDF-Projekt. Ich habe das damals ähnlich gemacht aber 
mit BPSK Kodierung gearbeitet. Du benutzt die MLS direkt. Durch die BPSK 
Kodierung wird zwar die Umlaufzeit der Sequenz verdoppelt aber das 
Signal auch im Frequenzbereich schmalbandiger. Da die BPSK 
Mittenfrequenz bei 24kHz lag wird so die Bandbreite des Soundsystems 
besser ausgenutzt. Zusätzlich wird das dekodierte Signal 
gleichspannungsfrei (Entfernung von tiefen Störfrequenzen eg. 50Hz 
Brumm) vor der FHT und das SNR verdoppelt.
Um nun das empfangene Signal wieder per FHT zu korrelieren muß man zwei 
FHTs um 90 Grad phasenverschoben durchführen. Man dekodiert also erstmal 
Sample für Sample die BPSK Kodierung und reduziert so auf die halbe 
Anzahl von Samples, zB. 2*511 Samples auf eine 511 Bit MLS. Man hat nach 
der Kodierung also zwei Blöcke a 511 Samples die mit 2 FHTs korreliert 
werden. Im Idealfall bekommt man bei der ersten FHT einen Peak mit 
511'er Amplitude und bei der zweiten FHT mit 0'er Amplitude. Wenn nun 
beide FHTs als Antwort 255 liefern als Peak dann kannst du diese 
Amplitudenwerte umrechnen in eine Phasenlage des Signal. Und das heist 
du hast so ein Verfahren das die Auflösung des Meßverfahren enorm 
erhöhen kann.

BPSK heist:
- ist das MLS Bit eine 0 so sende als Samples -1 und +1
- ist das MLS Bit eine 1 so sende als Samples +1 und -1

Bei der Dekodierung empfängst du doppelt so viele Samples wie Bits in 
der MLS. Subtrahiere nun jedes gerade Samples von nachfolgenden 
ungeradem Sample. Aus zb. -x und +x wird 2*x, ergo +3dB SNR. Aus diesem 
Datenstrom die geraden und ungerade Samples als zwei Ströme extrahieren 
und mit zwei FHTs korrelieren. Danach höchsten Peak in beiden Resultaten 
suchen und deren Peakwert + Index des Peaks in Relation Wellenlänge 
verrechnen wie oben schon angesprochen und schon hast du die Distanz. 
Beide Distanzen von Linken und rechten Kanal nur noch geometrisch zur 
Distanz der Lautsprecher verechnen und fertig ;)

Gruß hagen

[Detlef _. (detlef_a)]

Hallo Hagen,

>>>>>
Wochen beschäftigt ? Im Attachment alle 7710 primitiven Polynome für
Degree 17. Die Berechnung hat bei mir 1,3 Sekunden gedauert ;) Alle
276480 Polynome mit zB. Degree 24 zu finden dauert dann schon 9,7
Minuten.
<<<<

Gerade nachgekuckt, war fertig, ich habe auch 7710 gefunden :-)), davon 
4 mit nur zwei Rückkopplungen. Ich muß alle möglichen Polynome 
ausprobieren, weil ichs nicht besser verstehe, ich arbeite aber dran.

>>als Peak dann kannst du diese Amplitudenwerte umrechnen in eine Phasenlage des 
Signal.

Damit kommst du auf eine Auflösung besser als ein sample, geschickt. 
Deine BPSK Demodulation verstehe ich, aber das "muß man zwei FHTs um 90 
Grad phasenverschoben durchführen" verstehe ich nicht. Du schiebst die 
Folgen der geraden und der ungeraden samples durch diesselbe FHT, 
richtig!?

Sehr interessant, probier ich aus.

THX
Cheers
Detlef

[Detlef _. (detlef_a)]

Nachtrag, noch was zu primitiven Polynomen und 'trinomials', i.e. 
Schiebergister mit nur zwei Anzapfungen. Das ganze hängt mit Mersenne 
Primzahlen zusammen.

Cheers
Detlef

[Harald M. (mare_crisium)]

hagen,
die Idee mit der BPSK ist sehr schick, jedenfalls viel eleganter als 
meine Holzhammer- Modulation.

@detlef,
nee, Hagen teilt die 2*L Original-Samples in zwei Gruppen von je L 
Samples auf: In der ersten landet die Differenz der geraden und 
ungeraden, in der zweiten (nehme ich an) die Summe. Anschliessend 
schiebt er beide Gruppen einzeln durch die komplette Mühle mit 
Eingangs-Permutation, FHT und Ausgangspermutation. Das lefert dann das, 
was bei der Fouriertransformation "Amplituden-" und "Phasen-Spektrum" 
ist.

Let the good math roll!

Ciao,
Harald

[Hagen R. (hagen)]

Detlef _a schrieb:
> Damit kommst du auf eine Auflösung besser als ein sample, geschickt.
> Deine BPSK Demodulation verstehe ich, aber das "muß man zwei FHTs um 90
> Grad phasenverschoben durchführen" verstehe ich nicht. Du schiebst die
> Folgen der geraden und der ungeraden samples durch diesselbe FHT,
> richtig!?

Mal ein Beispiel:

MLS Sequenz lautet: 1 0 1 1 0 0 1
BPSK ergibt sich zu: +1 -1, -1 +1, +1 -1, +1 -1, -1 +1, -1 +1, +1 -1

Also doppelt so viele Samples wie Bits der MLS. Und der maximale Anzahl 
nachfolgender +1 und -1 Samples ist 2 der minimale ist 1. Wenn du das 
über die Soundkarte mit 24Khz sendest so haben wir Frequenzanteile von 
12kHz und 6Khz. Somit also viel schmalbandiger als wenn man die MLS 
direkt senden würde und d.h. besser angepasst an den Übertragungsbereich 
des Soundsystems.

Bei Empfang dieser BPSK Sequenz wird jedes Sample mit dem nachfolgenden 
subtrahiert und ergibt einen neuen Datenstrom.

Aus +1 -1 -1 +1 +1 -1 +1 -1 -1 +1 -1 +1 +1 -1 wird:

+2 0 -2 0 +2 -2 +2 0 -2 +2 -2 0 +2 -2

Nun diesen Strom in zwei Datenströme teilen, nach geredem und ungeradem 
Sample ergibt:

+2 -2 +2 +2 -2 -2 +2 = 1 0 1 1 0 0 1 unsere MLS

und

0  0 -2  0 +2  0 -2

diese Werte per FHT korrelieren und wir sehen das der obere Datenstrom 
unser MLS zu 100% enthält und der untere Datenstrom nicht.

Wenn dies zutrifft dann bist du Phasengleich zur MLS.

Normalerweise empfangen wir aber nicht +1 -1 Samples sondern +x und -x 
Samples mit unterschiedlicher Amplitude und auch DC Offset. Durch die 
Subtraktion zur Entfernung der BPSK Kodierung wird der DC-Offset 
eleminiert und damit auch Frequenzen ausserhalb unserer MLS Sequenz zb: 
50Hz Brumm. Ein DC Offset stellt für die Auswertung der FHT ein Problem 
dar bei der Ermittlung des Index des höchsten Peaks.

90 Grad ? weil wir 1 MLS Bit in zwei Samples aufteilen -> +1 -1, die 
Maximalamplituden einer Sinuswelle. Nach der Dekodierung der BPSK teilen 
wir den Datenstrom in zwei Teile. Im obigen Beispiel sind wir mit 0 Grad 
Phasenlage synchron beim ersten Datenstrom beim zweiten aber mit 90 Grad 
daneben und somit gibt es in diesem Datenstrom keinen Peak. Wären wir 
180 Grad Phasenversetzt so würde der Peak ins negative ausschlagen, wir 
hätten also eine invertierte MLS dekodiert.

Angenommen wird hätten ein 90 Grad Phasenverschobenes Signal, dann 
würden wir bei beiden Datenströmen einen Peak ins positive nach der FHT 
ermitteln der bei beiden zB. +3 und +4 betragen würde, summiert +7 wie 
bei einer 7 Bit MLS zu erwarten. Wäre er aber 270 Grad verschoben dann 
würden die Peaks bei -4 und -3 liegen. Mit Magnitude = 
Wurzel(Peak-Amplitude-1**2 + Peak-Amplitude-2**2) ermittelst du den 
Amplitudenwert und mit ArcTan() wie gewohnt die Phasenlage. Damit kannst 
du den ermittelten Indize der Peaks, der ja ganzzahlig ist in einen 
Realwert umrechnen und hast die Auflösung der Distanz-/Positionsmessung 
verbessert.

Gruß hagen

[Hagen R. (hagen)]

nun vergleichen wir mal beide Verfahren.

SNR ist gleich wenn: BPSK-MLS Lange = Direkt-MLS Länge / 2 ist. D.h. um 
das gleiche SNR zu bekommen benutzen wir einmal mit BPSK eine zB. 2^12-1 
lange MLS und bekommen eine BPSK Sequenz von 2^13-2. Die Verdoppeltung 
einer MLS Sequenz entspricht einer Anhebung des SNRs um +3dB also 2x 
besser. Wir nutzen also eine MLS direkt die eine Länge von 2^13-1 Bits 
hat damit sie vergleichbar zur BPSK-MLS wird. Das SNR ist also dann 
gleich bei beiden Verfahren und auch die Umlaufzeit der Sequenzen 
entspricht sich fast, +1 Samples weniger/mehr.

Aber: die BPSK-Dekodierung ist schneller als die FHT und die zwei FHTs 
sind nur marginal langsammer als eine FHT mit doppelter Größe. Auch die 
Zahlengrenzen der kürzeren FHT können noch im 16 Bit Bereich liegen 
während die doppelt so große FHT schon Überläufe erzeugen könnte (1 
Bit). Nutzt man diese Tatsache aus so heist dies das bei richtiger Wahl 
der MLS Länge die eine FHT die mit 16Bit Integern arbeitet viel 
performanter ist als die FHT die dann schon mit 17Bit Integern ergo 32 
Bit Integern arbeiten muß um den gleichen SNR erreichen zu können. Die 
Performancetechnische Umsetzung mit der BPSK Kodierung bringt also bei 
korrekter Wahl ein besseres Ergebnis. Gleiches gilt für den 
Speicherbedarf der gleichmal um das doppelte ansteigen würde.
Verbleibt noch der Vorteil das wir den DC-Offset entfernen und mit der 
IQ Demodulaton das Endergebnis noch höher auflösen können, für die 
BPSK-MLS Kodierung.

In Summa also ideal für kleine MCUs.

Gruß Hagen