3 van 4

Re: 10 Hz in het oor

Geplaatst: ma 02 sep 2024, 11:38
door Xilvo
wnvl1 schreef: ma 02 sep 2024, 11:25
Xilvo schreef: zo 01 sep 2024, 21:40 Dit krijg je (f0=10 Hz, samengesteld t.m. de duizendste harmonische, links zonder, rechte met willekeurige fase.
blok10Hz.png
De vorm blijft hier wel redelijk goed behouden.
Meer dan redelijk, sterker, die is perfect periodiek. ;)
Dat is logisch, ik transformeer niet maar genereer het signaal direct door de verschillende harmonische bij elkaar op te tellen. De faseverschuiving ligt per harmonische vast, dan moet het signaal wel periodiek worden.

Re: 10 Hz in het oor

Geplaatst: ma 02 sep 2024, 11:43
door wnvl1
1024, opnieuw niet perfect, maar wel redelijk.
Ook veel meer vervormd dan jullie figuren.
De reconstructie werkt perfect als ik niets verander aan de fase.
De fase veranderen doe ik door te vermenigvuldigen met een willekeurig punt op de eenheidscirkel in het complexe vlak voor de eerste helft en dan te vermenigvuldigen met het toegevoegde punt in de tweede helft.

Code: Selecteer alles

import numpy as np
import scipy 
from scipy import signal
import matplotlib.pyplot as plt

n=1024
t = np.linspace(0, 1, n)

zaagtand = signal.square(2 * np.pi * 5 * t)
plt.plot(t, zaagtand)
plt.title("blokgolf")

zaagtand_fft = scipy.fft.fft(signal.sawtooth(2 * np.pi * 5 * t)) 
random_fases = np.exp(1j *np.random.rand(n // 2 + 1) * 2 * np.pi)  
random_fases_gespiegeld = np.concatenate([random_fases, np.conj(random_fases[-2:0:-1])])
zaagtand_fft_gewijzigde_fase = zaagtand_fft * random_fases_gespiegeld
zaagtand_reconstructie = scipy.fft.ifft(zaagtand_fft_gewijzigde_fase)
plt.plot(t, zaagtand_reconstructie)
plt.title("blokgolf willekeurige fase")
download3

Re: 10 Hz in het oor

Geplaatst: ma 02 sep 2024, 11:45
door wnvl1
Xilvo schreef: ma 02 sep 2024, 11:38 Meer dan redelijk, sterker, die is perfect periodiek. ;)
Dat is logisch, ik transformeer niet maar genereer het signaal direct door de verschillende harmonische bij elkaar op te tellen. De faseverschuiving ligt per harmonische vast, dan moet het signaal wel periodiek worden.
Mijn opmerking was anders bedoeld. Ik had het niet over de periodiciteit, maar over de gelijkenis in vorm tussen één blok en wat jij bekomt na transformatie van de fases.

Re: 10 Hz in het oor

Geplaatst: ma 02 sep 2024, 11:51
door Xilvo
wnvl1 schreef: ma 02 sep 2024, 11:45 Mijn opmerking was anders bedoeld. Ik had het niet over de periodiciteit, maar over de gelijkenis in vorm tussen één blok en wat jij bekomt na transformatie van de fases.
Dat is toeval, dat hangt voornamelijk af van welke faseverschuiving de lagere harmonischen (derde, vijfde) krijgen.
Een andere keer, met andere willekeurige faseverschuivingen, ziet het er zo uit:
blok10Hz_a

Re: 10 Hz in het oor

Geplaatst: ma 02 sep 2024, 11:53
door HansH
wnvl1 schreef: ma 02 sep 2024, 11:22 Zelfde code, maar nu met een sampling van 1000 punten ipv 100punten.
bedoel je van 100 naar 1000 samples over de periode van 1s?
met 100 samples per 1 seconde is de sample afstand 10ms dus dan kun je frequentie componenten van 50 Hz of meer al niet meer weergeven dus als je die optelt krijg je fouten.
wat FFT daarmee te maken heeft zie ik nog even niet want je bent voor de signaal reconstructie sinussen aan het optellen. met FFT bereken je frequentiecomponenten, maar dat is hier niet nodig want je weet immers al uit de theorie welke componenten dat zijn en welke amplitude en fase.

Re: 10 Hz in het oor

Geplaatst: ma 02 sep 2024, 12:15
door wnvl1
HansH schreef: ma 02 sep 2024, 11:53 bedoel je van 100 naar 1000 samples over de periode van 1s?
ja
HansH schreef: ma 02 sep 2024, 11:53 met 100 samples per 1 seconde is de sample afstand 10ms dus dan kun je frequentie componenten van 50 Hz of meer al niet meer weergeven dus als je die optelt krijg je fouten.
maar de periodiciteit zou behouden moeten blijven
HansH schreef: ma 02 sep 2024, 11:53 wat FFT daarmee te maken heeft zie ik nog even niet want je bent voor de signaal reconstructie sinussen aan het optellen. met FFT bereken je frequentiecomponenten, maar dat is hier niet nodig want je weet immers al uit de theorie welke componenten dat zijn en welke amplitude en fase.
FFT is een 'luie' oplossing. Als je FFt gebruikt hoef je zelf niet meer de integralen uit te rekenen om de Fourier reeks op te stellen. In scipy zijn alle klassieke functies uit de signaalverwerking ingebouwd. Overigens wat een FFT doet is in de grond de Fourierreeks uitrekenen.

Re: 10 Hz in het oor

Geplaatst: ma 02 sep 2024, 12:21
door HansH
wnvl1 schreef: ma 02 sep 2024, 12:15
maar de periodiciteit zou behouden moeten blijven
dat klopt. je zou eens kunnen kijken als je alleen 1 of enkele van de frequentie componenten bekijkt als deel van de sommatie. dan moet je kunnen zien waar het probleem ontstaat.

Re: 10 Hz in het oor

Geplaatst: ma 02 sep 2024, 12:24
door HansH
wnvl1 schreef: ma 02 sep 2024, 12:15
FFT is een 'luie' oplossing. Als je FFt gebruikt hoef je zelf niet meer de integralen uit te rekenen om de Fourier reeks op te stellen. In scipy zijn alle klassieke functies uit de signaalverwerking ingebouwd. Overigens wat een FFT doet is in de grond de Fourierreeks uitrekenen.
dat is duidelijk, maar ik bedoelde dat de uitkomt van die integralen al gedaan is voor jou en je zo de formule voor elke frequentiecomponent beschikbaar hebt van de standaard vormen zoals blok, zaag en driehoek. (zoals xilvo ook al aangaf) zo heb ik het ook in LTspice grstopt uitgaande van de formule voor een blok

Re: 10 Hz in het oor

Geplaatst: ma 02 sep 2024, 12:57
door wnvl1
Mijn insteek is dat ik graag dit soort van problemen oplos met als uitgangspunten (1) zo snel mogelijk/zo weinig mogelijk lijnen code, (2) zo flexibel/algemeen mogelijk, (3) gebruik maken van alle libraries die al beschikbaar zijn en (4) als ik leuke wiskunde kan toepassen, kan dat het vorige zeker overrulen.

Bij deze opgave gingen mijn gedachte daarom onmiddellijk uit naar FFT.

Re: 10 Hz in het oor

Geplaatst: ma 02 sep 2024, 13:07
door HansH
dus je gebruikt de FFT om de frequentiecomponenten te berekenen en daarna tel je die componenten als functie van de tijd (of sample momenten in jouw geval) op om weer het oorspronkelijke signaal te resonstrueren als ik het goed begrijp. en als alternatief dezelfde frequentie componenten die uit de FFT rollen, maar dan met een random fase voor elke component?

Re: 10 Hz in het oor

Geplaatst: ma 02 sep 2024, 13:24
door wnvl1
Dat klopt. Dit zijn de vijf lijnen code waarin dit gebeurt.

Code: Selecteer alles

zaagtand_fft = scipy.fft.fft(signal.sawtooth(2 * np.pi * 5 * t)) 
random_fases = np.exp(1j *np.random.rand(n // 2 + 1) * 2 * np.pi)  
random_fases_gespiegeld = np.concatenate([random_fases, np.conj(random_fases[-2:0:-1])])
zaagtand_fft_gewijzigde_fase = zaagtand_fft * random_fases_gespiegeld
zaagtand_reconstructie = scipy.fft.ifft(zaagtand_fft_gewijzigde_fase)

Re: 10 Hz in het oor

Geplaatst: ma 02 sep 2024, 13:27
door Xilvo
De IFFT geeft complexe waardes, die plot je vervolgens. Gaat dat wel goed?

Re: 10 Hz in het oor

Geplaatst: ma 02 sep 2024, 13:34
door HansH
Xilvo schreef: ma 02 sep 2024, 13:27 De IFFT geeft complexe waardes, die plot je vervolgens. Gaat dat wel goed?
dat zijn als het goed is 2 toegevoegd complexe waardes die samen een reele amplitude en fase geven.
hier een voorbeeldje met fourier.
fourierheentrug
(455.25 KiB) 28 keer gedownload

Re: 10 Hz in het oor

Geplaatst: ma 02 sep 2024, 13:37
door wnvl1
Het complexe deel is heel klein, dat zou moeten vallen in de categorie numerieke fouten. Ik heb bvb na IFFT

[0.64439172-0.00197186j 0.76540039-0.0005845j 0.55170809-0.00197186j ...
0.69290493-0.0005845j 0.75350153-0.00197186j 0.4890989 -0.0005845j ]

Door van de random fase factoren in het tweede deel van de array die uit de FFT komt het complex toegevoegde te nemen, behoud je 'in principe' een reëel signaal. Dat gaat natuurlijk nooit helemaal exact kloppen.

Re: 10 Hz in het oor

Geplaatst: ma 02 sep 2024, 14:22
door Xilvo
wnvl1 schreef: ma 02 sep 2024, 13:37 Het complexe deel is heel klein, dat zou moeten vallen in de categorie numerieke fouten. Ik heb bvb na IFFT

[0.64439172-0.00197186j 0.76540039-0.0005845j 0.55170809-0.00197186j ...
0.69290493-0.0005845j 0.75350153-0.00197186j 0.4890989 -0.0005845j ]

Door van de random fase factoren in het tweede deel van de array die uit de FFT komt het complex toegevoegde te nemen, behoud je 'in principe' een reëel signaal. Dat gaat natuurlijk nooit helemaal exact kloppen.
Maar dit klopt wel erg slecht, dat zijn geen afrondingsfouten meer.

Component nul van de transform is het gemiddelde, daar is fase betekenisloos.
Component n/2 moet zijn eigen geconjugeerde zijn, moet reëel zijn.
Als ik beide componenten in random_fases_gespiegeld 1 maak, dan krijg ik in ieder geval een IFFT waarbij de imaginaire delen allemaal echt klein zijn, meestal kleiner dan 1E-15.