Forum » Programiranje » [Matlab] Fouriereva transformacija in analiza odzivov procesa
[Matlab] Fouriereva transformacija in analiza odzivov procesa
nosk8fx ::
Živjo,
v Matlabu moram narediti Fourierevo analizo odzivov realnega procesa. Se pravi, Fourierevo transformacijo vhoda (u) in izhoda (y), iz tega frekvenčni odziv, amplitudni in fazni del.
Od procesa imam tudi prenosno funkcijo (katera se zelo dobro prilega realnemu odzivu) in frekvenčni odziv te prenosne funkcije bi moral biti podoben frekvenčnemu odzivu realnih podatkov.
Najprej lažji del ...
frekvenčni odziv prenosne funkcije dobim s funkcijo bode:
[ampl_p,faza_p]=bode(syspd);%bode prenosne
amp_p=ampl_p(1,:);%amplitudni
faz_p=faza_p(1,:);%fazni
ampdb_p=20*log10(amp_p);%amplitudni v dB
nato pa max frekvenco, število vzorcev in vektor frekvenc:
fmax=1/2/ts;%max frekvenca
n_p=length(faz_p);%število vzorcev
frek_p=(0:1/ts/n_p/2:(fmax-1/ts/n_p/2))';%potek frekvence
pol to narišem in dobim v redu grafe.
Nadalje težji del ...
čas vzorčenja vzamem 50 s, zato iz realnih odzivov poberem vsak 50. vzorec, nato pa določim število vzorcev tukaj, resonančno frekvenco, max frekvenco in vektor frekvenc:
n=length(y);%število vzorcev
fres=1/n/ts;%resonančna
fmax=1/2/ts;%max
frek=(0:1/ts/n:fmax-1/ts/n)';%potek
potem naredim hitro Fourierevo transformacijo obeh in odstranim prvi vzorec:
ftu=fft(u,n);%vhoda
ftu(1)=[];
fty=fft(y,n);%identificirane
fty(1)=[];
nato po elementih zdelim obe transformaciji in vzamem samo prvo polovico signala, ker sta zrcalna:
g=fty./ftu;%frekvenčni odziv
for i=1:(n/2)%samo polovica signala
Ga(i,1)=g(i,1);
end
nazadnje še samo izračunam amplitudni (v dB) in fazni spekter (v stopinjah):
ampdb=20*log10(abs(Ga));
faz=180*unwrap(angle(Ga))/pi;
Vendar ko to narišem, so slike zelo čudne. Probaval sem tudi z različnimi časi vzorčenja in slike se seveda spreminjajo, vendar ni pravo.
Napaka mora biti nekje v tem algoritmu za Fourierevo transformacijo - mogoče sem kakšno pomembno podrobnost pozabil ali so mogoče vektorji (frekvenc ali česa drugega) narobe sestavljeni, vendar sam več nevem kaj naj.
Če se je mogoče kdo ukvarjal s tem oz. to dobro pozna, bi ga prosil za pomoč.
Upam, da nisem česar pozabil napisati.
Za vse že v naprej hvala.
Lep pozdrav,
Sašo
v Matlabu moram narediti Fourierevo analizo odzivov realnega procesa. Se pravi, Fourierevo transformacijo vhoda (u) in izhoda (y), iz tega frekvenčni odziv, amplitudni in fazni del.
Od procesa imam tudi prenosno funkcijo (katera se zelo dobro prilega realnemu odzivu) in frekvenčni odziv te prenosne funkcije bi moral biti podoben frekvenčnemu odzivu realnih podatkov.
Najprej lažji del ...
frekvenčni odziv prenosne funkcije dobim s funkcijo bode:
[ampl_p,faza_p]=bode(syspd);%bode prenosne
amp_p=ampl_p(1,:);%amplitudni
faz_p=faza_p(1,:);%fazni
ampdb_p=20*log10(amp_p);%amplitudni v dB
nato pa max frekvenco, število vzorcev in vektor frekvenc:
fmax=1/2/ts;%max frekvenca
n_p=length(faz_p);%število vzorcev
frek_p=(0:1/ts/n_p/2:(fmax-1/ts/n_p/2))';%potek frekvence
pol to narišem in dobim v redu grafe.
Bodejev diagram prenosne funkcije
Nadalje težji del ...
čas vzorčenja vzamem 50 s, zato iz realnih odzivov poberem vsak 50. vzorec, nato pa določim število vzorcev tukaj, resonančno frekvenco, max frekvenco in vektor frekvenc:
n=length(y);%število vzorcev
fres=1/n/ts;%resonančna
fmax=1/2/ts;%max
frek=(0:1/ts/n:fmax-1/ts/n)';%potek
potem naredim hitro Fourierevo transformacijo obeh in odstranim prvi vzorec:
ftu=fft(u,n);%vhoda
ftu(1)=[];
fty=fft(y,n);%identificirane
fty(1)=[];
nato po elementih zdelim obe transformaciji in vzamem samo prvo polovico signala, ker sta zrcalna:
g=fty./ftu;%frekvenčni odziv
for i=1:(n/2)%samo polovica signala
Ga(i,1)=g(i,1);
end
nazadnje še samo izračunam amplitudni (v dB) in fazni spekter (v stopinjah):
ampdb=20*log10(abs(Ga));
faz=180*unwrap(angle(Ga))/pi;
Vendar ko to narišem, so slike zelo čudne. Probaval sem tudi z različnimi časi vzorčenja in slike se seveda spreminjajo, vendar ni pravo.
Napaka mora biti nekje v tem algoritmu za Fourierevo transformacijo - mogoče sem kakšno pomembno podrobnost pozabil ali so mogoče vektorji (frekvenc ali česa drugega) narobe sestavljeni, vendar sam več nevem kaj naj.
Če se je mogoče kdo ukvarjal s tem oz. to dobro pozna, bi ga prosil za pomoč.
Upam, da nisem česar pozabil napisati.
Za vse že v naprej hvala.
Lep pozdrav,
Sašo
rasta ::
Najprej razčistimo izrazoslovje. "max frekvenca" naj bi bila frekvenca zgibanja, "resonančna frekvenca" pa prirastek frekvence?
Nadalje, "v redu grafi" meni izgledajo čudni, sploh ker so analitično določeni (krivulja bi morala biti gladka, brez artefaktov). Je prenosna funkcija linearna in katerega reda je (verjetno drugega)? Prenosno funkcijo si pa določil z identifikacijo z najmanjšimi kvadrati?
Koliko je število vzorcev n? Je potenca števila 2?
Če ni, hitra fourierova transformacija doda toliko ničel, da je. Problem pri tem pa je predvsem, ko signal ob končnem času ni približno nič, saj takrat dobiš oster rob in veliko artifaktov v frekvenčnem spektru.
Lahko malo bolje opišeš, na kakšen način so "čudne".
Kakšna pa je Nyquistova frekvenca zvorčenega procesa?
Zakaj po izvedbi hitre Fourierove transformacije odvzameš prvi vzorec (statično ojačenje)?
Lahko mogoče navedeš še, kako narišeš slike?
Nadalje, "v redu grafi" meni izgledajo čudni, sploh ker so analitično določeni (krivulja bi morala biti gladka, brez artefaktov). Je prenosna funkcija linearna in katerega reda je (verjetno drugega)? Prenosno funkcijo si pa določil z identifikacijo z najmanjšimi kvadrati?
Koliko je število vzorcev n? Je potenca števila 2?
Če ni, hitra fourierova transformacija doda toliko ničel, da je. Problem pri tem pa je predvsem, ko signal ob končnem času ni približno nič, saj takrat dobiš oster rob in veliko artifaktov v frekvenčnem spektru.
"Vendar ko to narišem, so slike zelo čudne."
Lahko malo bolje opišeš, na kakšen način so "čudne".
Kakšna pa je Nyquistova frekvenca zvorčenega procesa?
Zakaj po izvedbi hitre Fourierove transformacije odvzameš prvi vzorec (statično ojačenje)?
Lahko mogoče navedeš še, kako narišeš slike?
nosk8fx ::
ja, izrazi so tako nekako "isti" 
...
glede grafov z bodeja, sem jih malo popravil, ker sem uporabil vektor frekvence od funkcije bode, prej pa sem vektor sestavil z linspace (vektor od funkcije bode ne gre z enakomernimi prirastki).
tudi algoritem FT sem malo dopolnil (po Matlab helpu od fft), najprej pa še vprašanje čisto od začetka:
signal z realnega procesa ima nekje 22000 vzorcev in je samplan na 1 s, vendar je to za identifikacijo preveč, pa tudi z metodo najmanjših kvadratov sem delal s ts=50 s in čisto lepo dela. se pravi, moram za FT pobrat vsak 50. vzorec in potem imam krajši y in u.
nato je število vzorcev enako dolžini y, pri ts=50 oz. fs=1/50:
n=length(y);%število mojih vzorcev
nato pa število vzorcev za FT povečam na naslednjo potenco 2, naredim FT (pustim vse vzorce) in jo normiram z številom mojih vzorcev:
N=2^nextpow2(n); % naslednja potenca 2 od vektorja y
fty=fft(y,N)/n;%FT vhoda
ftu=fft(u,N)/n;%FT izhoda
za vektor frekvence pa naredim linspace pomnožen s fs/2:
frek=1/ts/2*linspace(0,1,N/2+1);%vektor frekvence
pol zdelim fty/ftu da dobim frekvenčno karakteristiko, izračunam amplitudni del (v dB) in fazni (v stopinjah) samo prve polovice signala (mogoče veš, zakaj v Matlab helpu ko riše amplitudni del, da 2*abs? - jaz sem zdaj risal brez 2):
g=fty./ftu;%frekvenčna karakteristika
ampdb=20*log10(abs(g(1:N/2+1)));%amplitudni
faz=unwrap(180/pi*angle(g(1:N/2+1)));%fazni
nato pa te grafe narišem z semilogx, vendar so slike čudne ...
najbolj čudni so grafi FT samo u-ja ... u je pri meni stopnica z 0 na 5500 ... res je, da to ni najboljši signal, samo drugega nimam.
graf y še dobro zgleda, samo u je pa katastrofa in pol je seveda tudi celotna karakteristika narobe.
kakšen predlog kaj je narobe pri u?
pa glede vektorja frekvence, ki ga sestavljam sam - v Matlab helpu je narjen od 0 z enakomernimi koraki, vendar pa je pol naslednji pri 10^-5, vektor od funkcije bode pa začne pri 10^-7 in v log merilu je to kar velika razlika. kakšen predlog kako to nekako "uskladit"?
hvala in lp,
Sašo
...glede grafov z bodeja, sem jih malo popravil, ker sem uporabil vektor frekvence od funkcije bode, prej pa sem vektor sestavil z linspace (vektor od funkcije bode ne gre z enakomernimi prirastki).
boljši grafi
tudi algoritem FT sem malo dopolnil (po Matlab helpu od fft), najprej pa še vprašanje čisto od začetka:
signal z realnega procesa ima nekje 22000 vzorcev in je samplan na 1 s, vendar je to za identifikacijo preveč, pa tudi z metodo najmanjših kvadratov sem delal s ts=50 s in čisto lepo dela. se pravi, moram za FT pobrat vsak 50. vzorec in potem imam krajši y in u.
nato je število vzorcev enako dolžini y, pri ts=50 oz. fs=1/50:
n=length(y);%število mojih vzorcev
nato pa število vzorcev za FT povečam na naslednjo potenco 2, naredim FT (pustim vse vzorce) in jo normiram z številom mojih vzorcev:
N=2^nextpow2(n); % naslednja potenca 2 od vektorja y
fty=fft(y,N)/n;%FT vhoda
ftu=fft(u,N)/n;%FT izhoda
za vektor frekvence pa naredim linspace pomnožen s fs/2:
frek=1/ts/2*linspace(0,1,N/2+1);%vektor frekvence
pol zdelim fty/ftu da dobim frekvenčno karakteristiko, izračunam amplitudni del (v dB) in fazni (v stopinjah) samo prve polovice signala (mogoče veš, zakaj v Matlab helpu ko riše amplitudni del, da 2*abs? - jaz sem zdaj risal brez 2):
g=fty./ftu;%frekvenčna karakteristika
ampdb=20*log10(abs(g(1:N/2+1)));%amplitudni
faz=unwrap(180/pi*angle(g(1:N/2+1)));%fazni
nato pa te grafe narišem z semilogx, vendar so slike čudne ...
frekvenčna karakteristika
samo FT y
samo FT u
najbolj čudni so grafi FT samo u-ja ... u je pri meni stopnica z 0 na 5500 ... res je, da to ni najboljši signal, samo drugega nimam.
graf y še dobro zgleda, samo u je pa katastrofa in pol je seveda tudi celotna karakteristika narobe.
kakšen predlog kaj je narobe pri u?
pa glede vektorja frekvence, ki ga sestavljam sam - v Matlab helpu je narjen od 0 z enakomernimi koraki, vendar pa je pol naslednji pri 10^-5, vektor od funkcije bode pa začne pri 10^-7 in v log merilu je to kar velika razlika. kakšen predlog kako to nekako "uskladit"?
hvala in lp,
Sašo
nosk8fx ::
zdaj sem probal še eno drugo stvar ... ker je u itak samo konstanta, niti ni potrebno delati FT, saj je amplitudni spekter 20*log10(u), faze pa itak ne doda nič.
tako da sem za amplitudni del samo delil fty z 20*log10(u), fazo pa upošteval samo od fty.
dobim boljšo sliko, vendar ni še prava ...
amplitudni spekter je previsoko, pa verjetno dela probleme tudi vektor frekvence, ki je različen.
kakšni predlogi sedaj, kako to rešit?
lp
tako da sem za amplitudni del samo delil fty z 20*log10(u), fazo pa upošteval samo od fty.
dobim boljšo sliko, vendar ni še prava ...
amplitudni spekter je previsoko, pa verjetno dela probleme tudi vektor frekvence, ki je različen.
kakšni predlogi sedaj, kako to rešit?
lp
nosk8fx ::
sem v zadnjem postu naredil napako pri tem deljenju, delil sem tako:
ampdb=20*log10(abs(fty(1:N/2+1))./u(1:N/2+1));
slika pa je
se pravi, je amplitudni spekter prenizko, premaknjene frekvence pa ostanejo ...
lp
ampdb=20*log10(abs(fty(1:N/2+1))./u(1:N/2+1));
slika pa je
se pravi, je amplitudni spekter prenizko, premaknjene frekvence pa ostanejo ...
lp
rasta ::
Še vedno nisi jasno razložil, kaj počneš.
Če prav razumem imaš vektor izmerkov - odziv sistema na stopnico. Iz teh izmerkov najprej z metodo najmanjših kvadratov identificiraš prenosno funkcijo, ki je drugega reda z dvemi ničlami.
Nato sistem identificiraš še s Fourierovo analizo in primerjaš z odzivom ocenjene prenosne funkcije. Odziva se ne skladata in to te muči.
Kar se tiče decimacije in spremembe časa vzorčenja - glede na prvi Bodejev diagram (kako si dobil ta diagram) imaš eno ničlo pri frekvenci, ki je višja od polovice frekvence vzorčenja 1/50s, se pravi s tem pokvariš signal (-> aliasing).
Tu si naredil padding signala - podaljšal si ga z ničlami do N. Prav enako stvar ti naredi že metoda fft (se pravi je ta N čisto odveč).
Ni pa nujno, da to hočeš - če je zadnji vzorec signala zelo različen od nič, boš dobil v transformaciji kar precej artefaktov.
To popravljaš nek koeficient pri fft-ju, za kaj točno že gre pa sem že pozabil.
Kateri graf pa prikazuje FFT u-ja? Pa zakaj je v logaritmskem merilu?
Bodejev diagram prenosne funkcije (katero imaš podano analitično) zračunaš pri custom frekvencah (poglej si help od funkcije bode), najbolj zaželjeno je seveda, da pri enakih kot pri neparametrični identifikaciji.
Nisi malo prej napisal, da je u stopnica? Takrat amplitudni spekter ni konstanta ...
Sploh pa če je u konstanta, boš imel kar precej težav z identifikacijo (se je v bistvu ne da izvest).
Tvoj opis problema in postopka je tako zmeden, da ne morem uganiti, kaj je narobe. Ne vem niti, kaj je na katerem od tvojih grafov. Če ni skrivnost, je verjetno najbolje, da objaviš kar celoten .m fajl.
Če prav razumem imaš vektor izmerkov - odziv sistema na stopnico. Iz teh izmerkov najprej z metodo najmanjših kvadratov identificiraš prenosno funkcijo, ki je drugega reda z dvemi ničlami.
Nato sistem identificiraš še s Fourierovo analizo in primerjaš z odzivom ocenjene prenosne funkcije. Odziva se ne skladata in to te muči.
Kar se tiče decimacije in spremembe časa vzorčenja - glede na prvi Bodejev diagram (kako si dobil ta diagram) imaš eno ničlo pri frekvenci, ki je višja od polovice frekvence vzorčenja 1/50s, se pravi s tem pokvariš signal (-> aliasing).
N=2^nextpow2(n); % naslednja potenca 2 od vektorja y
fty=fft(y,N)/n;%FT vhoda
ftu=fft(u,N)/n;%FT izhoda
Tu si naredil padding signala - podaljšal si ga z ničlami do N. Prav enako stvar ti naredi že metoda fft (se pravi je ta N čisto odveč).
Ni pa nujno, da to hočeš - če je zadnji vzorec signala zelo različen od nič, boš dobil v transformaciji kar precej artefaktov.
(mogoče veš, zakaj v Matlab helpu ko riše amplitudni del, da 2*abs?
To popravljaš nek koeficient pri fft-ju, za kaj točno že gre pa sem že pozabil.
kakšen predlog kaj je narobe pri u?
Kateri graf pa prikazuje FFT u-ja? Pa zakaj je v logaritmskem merilu?
kakšen predlog kako to nekako "uskladit"?
Bodejev diagram prenosne funkcije (katero imaš podano analitično) zračunaš pri custom frekvencah (poglej si help od funkcije bode), najbolj zaželjeno je seveda, da pri enakih kot pri neparametrični identifikaciji.
zdaj sem probal še eno drugo stvar ... ker je u itak samo konstanta, niti ni potrebno delati FT, saj je amplitudni spekter 20*log10(u), faze pa itak ne doda nič.
Nisi malo prej napisal, da je u stopnica? Takrat amplitudni spekter ni konstanta ...
Sploh pa če je u konstanta, boš imel kar precej težav z identifikacijo (se je v bistvu ne da izvest).
Tvoj opis problema in postopka je tako zmeden, da ne morem uganiti, kaj je narobe. Ne vem niti, kaj je na katerem od tvojih grafov. Če ni skrivnost, je verjetno najbolje, da objaviš kar celoten .m fajl.
nosk8fx ::
po vrsti ...
da, imam vektor izmerkov, ki je odziv sistema na stopnico. in ja, z metodo najmanjših kvadratov identificirana prenosna funkcija je drugega reda z dvema ničlama. nato pa jo moram primerjati še z odzivi, dobljenimi s F. analizo, in ti odzivi se ne skladajo, kot si napisal.
tole z ničlami in frekvenco vzorčenja sem malce spregledal, saj pride polovica 0,01 /s, ena ničla pa je (s+0,022) (druga je (s+0,00015)). se pravi moram čas vzorčenja spustit na nekje 22 oz. 20 s.
tole z velikim N je narejeno enako v Matlab helpu, zato sem še to poizkusil (na začetku sem delal brez). kako pol, število vzorcev n pol vzamem samo število dolžine vektorja y, se pravi n=length(y);, pol pa samo fty=fft(y,n); ali moram tudi normirati fty=fft(y,n)/n;?
se pravi, da ko računam absolutno vrednost, moram pomnožiti z 2 zaradi popravljanja koeficientov fft.
graf FFT u-ja prikazuje ta spodnji (zgoraj piše na desni strani slike):
v logaritemskem merilu sem risal, ker funkcija bode riše v logaritemskem na abscisi (rad/s), amplitudo v dB, fazo pa v stopinjah.
se pravi, da bi bilo boljše, da tudi pri funkciji bode uporabim isti vektor frekvenc (v defaultu označen kot w).
pa če u gledam kot stopnico, ki nastopi v času 0 (od 0 dalje gledam signale), je to konstanta - tako imam tudi sestavljen u, da so same iste cifre. vem, da je sam vhodni signal slab za identifikacijo, vendar bi se menda tukaj dalo tako, da ne delam FFT u-ja.
zdaj mislim, da bi moral razumeti, kaj me muči. prosim, poglej to, kar sem zdaj napisal, če je še kaj narobe. sam pa bom probal stvar identificirati še pri manjšem času vzorčenja, pa bom pol videl kako bo ratalo, upošteval pa bom tudi moj vektor frekvenc pri funkciji bode.
pol če bo treba, bom prilepil .m file, ni problema.
lp
da, imam vektor izmerkov, ki je odziv sistema na stopnico. in ja, z metodo najmanjših kvadratov identificirana prenosna funkcija je drugega reda z dvema ničlama. nato pa jo moram primerjati še z odzivi, dobljenimi s F. analizo, in ti odzivi se ne skladajo, kot si napisal.
tole z ničlami in frekvenco vzorčenja sem malce spregledal, saj pride polovica 0,01 /s, ena ničla pa je (s+0,022) (druga je (s+0,00015)). se pravi moram čas vzorčenja spustit na nekje 22 oz. 20 s.
tole z velikim N je narejeno enako v Matlab helpu, zato sem še to poizkusil (na začetku sem delal brez). kako pol, število vzorcev n pol vzamem samo število dolžine vektorja y, se pravi n=length(y);, pol pa samo fty=fft(y,n); ali moram tudi normirati fty=fft(y,n)/n;?
se pravi, da ko računam absolutno vrednost, moram pomnožiti z 2 zaradi popravljanja koeficientov fft.
graf FFT u-ja prikazuje ta spodnji (zgoraj piše na desni strani slike):
v logaritemskem merilu sem risal, ker funkcija bode riše v logaritemskem na abscisi (rad/s), amplitudo v dB, fazo pa v stopinjah.
se pravi, da bi bilo boljše, da tudi pri funkciji bode uporabim isti vektor frekvenc (v defaultu označen kot w).
pa če u gledam kot stopnico, ki nastopi v času 0 (od 0 dalje gledam signale), je to konstanta - tako imam tudi sestavljen u, da so same iste cifre. vem, da je sam vhodni signal slab za identifikacijo, vendar bi se menda tukaj dalo tako, da ne delam FFT u-ja.
zdaj mislim, da bi moral razumeti, kaj me muči. prosim, poglej to, kar sem zdaj napisal, če je še kaj narobe. sam pa bom probal stvar identificirati še pri manjšem času vzorčenja, pa bom pol videl kako bo ratalo, upošteval pa bom tudi moj vektor frekvenc pri funkciji bode.
pol če bo treba, bom prilepil .m file, ni problema.
lp
rasta ::
tole z velikim N je narejeno enako v Matlab helpu, zato sem še to poizkusil (na začetku sem delal brez).
Sem šel pogledat Matlabovo dokumentacijo za fft in ugotovil, da formulo za FFT uporabi samo takrat, ko je dolžina signala potenca števila dva (kar je pogoj), drugače pa uporabi neko drugo formulo (nekakšno pseudo-FFT).
Zanimivo. Človek se vsak dan kaj novega nauči.
se pravi, da ko računam absolutno vrednost, moram pomnožiti z 2 zaradi popravljanja koeficientov fft.
V bistu, če nisi opazil, rezultat FFT-ja deliš z L/2.
Ker sem FFT študiral (in kodiral) že kar nekaj let nazaj, se res ne spomnim več, kako je s temi faktorji.
Uporabiš pač pragmatičen, inženirski princip: pri nizkih frekvenca mora biti amplituda približno enaka statičnem ojačanju, če ni, pomnožiš s kakšnim faktorjem
pa če u gledam kot stopnico, ki nastopi v času 0 (od 0 dalje gledam signale), je to konstanta - tako imam tudi sestavljen u, da so same iste cifre.
Ustavimo se najprej pri vhodnem signalu.
Vemo, da je Fourierov transform stopnice delta impulz oz. pri diskretni Fourierovi transformaciji so to dolgi pravokotni pulzi (ker po definiciji transformira periodične signale). Fourier pravokotnih signalov je funkcija oblike sin(x)/x.
Iz Bodejevega grafa vhodnega signala je razvidno, da je prva nula že pri 1.4*10-4, do tu se tretira vzbujanje spodobno, pri višjih frekvencah pa boš za oceno dobil same smeti. Iz Bodejevega diagrama identificirane prenosne funkcije je razvidno, da je to ravno nekako mejna frekvenca sistema. Kaj dosti se torej na Bodejevem diagramu dobljenim z Fourierovo analizo ne bo videlo.
rasta ::
Pa še opozoril bi te, da Matlabova funkcija bode dela s krožnimi frekvencami (frekvenca*2*π) -- poglej, če pri plotanju slučajno kje ne mešaš ti dve frekvenci. (Ta "bug" dostikrat doprinese k zanimivim rezultatom ...)
nosk8fx ::
se pravi, postopek izračuna FFT je takole v redu:
n=length(y);
N=2^nextpow2(n);
fty=fft(y,N)/(n/2);
ja, res je, da je FT stopnice delta impulz oz. diskretno sin(x)/x. samo, če pa gledam Bodejev diagram (na katerem moram oba odziva primerjati) konstante, je to pri amplitudi tudi konstanta (20*log(K)), faze pa ne daje. zato sem poizkusil tudi direktno s tem deljenjem, pa mi daje "lepše" rezultate.
tole, ko si napisal glede krožnih frekvenc funkcije bode - tole nisem razumel najbolje. bi moral kje pomnožiti ali deliti ali kako da nebi mešal teh frekvenc?
no, sedaj sem poizkusil s časom vzorčenja 20 s, da je višja ničla nižja od polovice frekvence vzorčenja (zdaj če bi gledal nižjo ničlo, bi to prišli zelo niziki časi vzorčenja - res je, da bo od tiste frekvence naprej slabo, ampak naj se vidi vsaj približno ujemanje do tam):
prvo, s fft obeh (u in y), pa nato z deljenjem po elementih g=fty./ftu;, nato pa še absolutni in fazni del (ampdb=20*log10(abs(g(1:N/2+1))); in faz=unwrap(180/pi*angle(g(1:N/2+1)));). pri risanju s funkcijo bode sem mu "vsilil" to mojo sestavljeno frekvenco (frek=1/ts/2*linspace(0,1,N/2+1);). tako dobim:
drugo, pa sem fft nardil samo y, nato pa pri absolutni vrednosti delil s konstanto (ampdb=20*log10(abs(fty(1:N/2+1))./u(1:N/2+1));). tako pa dobim:
kako to komentirat, nevem točno ...
pri prvi sliki je amplitudni del ful narobe, čeprav sta sami FFT ok (sem pogledal grafe), pri drugi sliki pa je oblika amplitudnega dokaj v redu, vendar je malce prenizko.
pri prvi sliki pol faza na koncu preveč narase, ampak lahko da tistega dela več sploh ne gledamo, pri drugi sliki pa pol preveč pade, ampak je spet to že v območju "smeti".
lp
n=length(y);
N=2^nextpow2(n);
fty=fft(y,N)/(n/2);
ja, res je, da je FT stopnice delta impulz oz. diskretno sin(x)/x. samo, če pa gledam Bodejev diagram (na katerem moram oba odziva primerjati) konstante, je to pri amplitudi tudi konstanta (20*log(K)), faze pa ne daje. zato sem poizkusil tudi direktno s tem deljenjem, pa mi daje "lepše" rezultate.
tole, ko si napisal glede krožnih frekvenc funkcije bode - tole nisem razumel najbolje. bi moral kje pomnožiti ali deliti ali kako da nebi mešal teh frekvenc?
no, sedaj sem poizkusil s časom vzorčenja 20 s, da je višja ničla nižja od polovice frekvence vzorčenja (zdaj če bi gledal nižjo ničlo, bi to prišli zelo niziki časi vzorčenja - res je, da bo od tiste frekvence naprej slabo, ampak naj se vidi vsaj približno ujemanje do tam):
prvo, s fft obeh (u in y), pa nato z deljenjem po elementih g=fty./ftu;, nato pa še absolutni in fazni del (ampdb=20*log10(abs(g(1:N/2+1))); in faz=unwrap(180/pi*angle(g(1:N/2+1)));). pri risanju s funkcijo bode sem mu "vsilil" to mojo sestavljeno frekvenco (frek=1/ts/2*linspace(0,1,N/2+1);). tako dobim:
drugo, pa sem fft nardil samo y, nato pa pri absolutni vrednosti delil s konstanto (ampdb=20*log10(abs(fty(1:N/2+1))./u(1:N/2+1));). tako pa dobim:
kako to komentirat, nevem točno ...
pri prvi sliki je amplitudni del ful narobe, čeprav sta sami FFT ok (sem pogledal grafe), pri drugi sliki pa je oblika amplitudnega dokaj v redu, vendar je malce prenizko.
pri prvi sliki pol faza na koncu preveč narase, ampak lahko da tistega dela več sploh ne gledamo, pri drugi sliki pa pol preveč pade, ampak je spet to že v območju "smeti".
lp
Vredno ogleda ...
| Tema | Ogledi | Zadnje sporočilo | |
|---|---|---|---|
| Tema | Ogledi | Zadnje sporočilo | |
| » | Bodejev diagramOddelek: Šola | 3082 (1928) | čuhalev |
| » | Arduino frekvenca zvokaOddelek: Elektrotehnika in elektronika | 981 (854) | sandmat |
| » | Fourierjeva vrsta, fourierjeva transformacijaOddelek: Šola | 6916 (3702) | marjan_h |
| » | spektralni analizatorOddelek: Elektrotehnika in elektronika | 3420 (2471) | elektro87 |
| » | merjenje obratovOddelek: Programiranje | 1692 (1484) | neoto |