MATLAB FFT अप खिलवाड़ सीमा और fftshift

Question

xaxis यह पहली बार मैं fft समारोह का उपयोग कर रहा है और मैं एक साधारण कोसाइन फंक्शन के आवृत्ति स्पेक्ट्रम साजिश कोशिश कर रहा हूँ है:

च = क्योंकि (2 * pi * 300 * टी)

नमूना दर 220500 है। मैं फ़ंक्शन f का एक सेकंड प्लॉट कर रहा हूं।

यहाँ मेरी प्रयास है:

time = 1; 
freq = 220500; 
t = 0 : 1/freq : 1 - 1/freq; 
N = length(t); 
df = freq/(N*time); 

F = fftshift(fft(cos(2*pi*300*t))/N); 
faxis = -N/2/time : df : (N/2-1)/time; 

plot(faxis, real(F)); 
grid on; 
xlim([-500, 500]); 

जब मैं 900Hz करने के लिए आवृत्ति को बढ़ाने क्यों मैं अजीब परिणाम मिलता है? इन विषम परिणामों को एक्स-अक्ष सीमाओं को 500 हर्ट्ज से 1000 हर्ट्ज तक बढ़ाकर तय किया जा सकता है। साथ ही, क्या यह सही दृष्टिकोण है? मैंने देखा कि कई अन्य लोगों ने fftshift(X) का उपयोग नहीं किया है (लेकिन मुझे लगता है कि उन्होंने केवल एक तरफा स्पेक्ट्रम विश्लेषण किया था)।

धन्यवाद।

Answer 1

यहां वादा के रूप में मेरी प्रतिक्रिया है।

आपके पहले या आपके प्रश्न "आप आवृत्ति को 900 हर्ट्ज तक बढ़ाते समय अजीब परिणाम प्राप्त करते हैं" से संबंधित है, जैसा कि @ कैस्टिलो द्वारा वर्णित मैटलैब की साजिश रेजेलिंग कार्यक्षमता से संबंधित है। जब आप एक्स-अक्ष की सीमा बदलते हैं, तो मैटलैब सहायक होने का प्रयास करेगा और वाई-अक्ष को पुन: सहेज लेगा। यदि चोटियों को आपकी निर्दिष्ट सीमा से बाहर झूठ बोलती है, तो मैटलैब प्रक्रिया में उत्पन्न छोटी संख्यात्मक त्रुटियों पर ज़ूम इन करेगा। यदि आप परेशान करते हैं तो आप इसे 'यलिम' कमांड के साथ उपाय कर सकते हैं।

हालांकि, आपका दूसरा, अधिक खुला प्रश्न "क्या यह सही दृष्टिकोण है?" एक गहरी चर्चा की आवश्यकता है। मुझे आपको यह बताने की अनुमति दें कि मैं कोसाइन लहर की साजिश के अपने लक्ष्य को प्राप्त करने के लिए एक और अधिक लचीला समाधान बनाने के बारे में कैसे जाऊंगा।

आप निम्न के साथ शुरू:

time = 1; 
freq = 220500; 

यह तुरंत मेरे सिर में एक अलार्म को जन्म देती है। बाकी पोस्ट को देखते हुए, आपको उप-केएचजेड रेंज में आवृत्तियों में दिलचस्पी दिखाई देती है। यदि ऐसा है, तो यह नमूना दर अत्यधिक है क्योंकि इस दर के लिए Nyquist सीमा (sr/2) 100 kHz से ऊपर है।मैं अनुमान लगा रहा हूं कि 22050 हर्ट्ज की सामान्य ऑडियो नमूना दर का उपयोग करना है (लेकिन मैं यहां गलत हो सकता हूं)?

किसी भी तरह से, आपका विश्लेषण अंत में संख्यात्मक रूप से ठीक काम करता है। हालांकि, आप यह समझने में मदद नहीं कर रहे हैं कि असली दुनिया स्थितियों में विश्लेषण के लिए एफएफटी का सबसे प्रभावी ढंग से उपयोग कैसे किया जा सकता है।

मुझे यह पोस्ट करने की अनुमति दें कि मैं यह कैसे करूं। निम्न स्क्रिप्ट लगभग आपकी स्क्रिप्ट क्या करती है, लेकिन कुछ संभावित खुलती है जिस पर हम निर्माण कर सकते हैं। ।

%// These are the user parameters 
durT = 1; 
fs = 22050; 
NFFT = durT*fs; 
sigFreq = 300; 

%//Calculate time axis 
dt = 1/fs; 
tAxis = 0:dt:(durT-dt); 

%//Calculate frequency axis 
df = fs/NFFT; 
fAxis = 0:df:(fs-df); 

%//Calculate time domain signal and convert to frequency domain 
x = cos( 2*pi*sigFreq*tAxis ); 
F = abs( fft(x, NFFT)/NFFT ); 

subplot(2,1,1); 
plot( fAxis, 2*F ) 
xlim([0 2*sigFreq]) 
title('single sided spectrum') 

subplot(2,1,2); 
plot( fAxis-fs/2, fftshift(F) ) 
xlim([-2*sigFreq 2*sigFreq]) 
title('whole fft-shifted spectrum') 

आप समय अक्ष की गणना करते हैं और समय धुरी की लंबाई से आपकी एफएफटी अंकों की गणना करते हैं। यह बहुत अजीब है। इस दृष्टिकोण के साथ समस्या यह है कि जब आप अपने इनपुट सिग्नल की अवधि बदलते हैं तो एफएफटी के आवृत्ति रिज़ॉल्यूशन में परिवर्तन होता है, क्योंकि एन आपके "टाइम" चर पर निर्भर होता है। Matlab fft कमांड एक एफएफटी आकार का उपयोग करेगा जो इनपुट सिग्नल के आकार से मेल खाता है।

मेरे उदाहरण में, मैं सीधे एनएफएफटी से आवृत्ति अक्ष की गणना करता हूं। यह उपर्युक्त उदाहरण के संदर्भ में कुछ हद तक अप्रासंगिक है, क्योंकि मैंने सिग्नल में नमूने की संख्या के बराबर एनएफएफटी सेट किया है। हालांकि, इस प्रारूप का उपयोग करने से आपकी सोच को कम करने में मदद मिलती है और यह मेरे अगले उदाहरण में बहुत महत्वपूर्ण हो जाता है।

** साइड नोट: आप अपने उदाहरण में वास्तविक (एफ) का उपयोग करते हैं। जब तक आपके पास एफएफटी परिणाम के वास्तविक हिस्से को निकालने का बहुत अच्छा कारण न हो, तब तक पेट (एफ) का उपयोग करके एफएफटी की परिमाण निकालने के लिए यह बहुत आम है। यह sqrt (वास्तविक (एफ) के बराबर है।^2 + कल्पना (एफ)।^2)। **

अधिकांश समय आप एक छोटे एनएफएफटी का उपयोग करना चाहते हैं। ऐसा इसलिए हो सकता है क्योंकि आप शायद वास्तविक समय प्रणाली में विश्लेषण चला रहे हैं, या क्योंकि आप एक समय के लिए अलग-अलग स्पेक्ट्रम का विचार प्राप्त करने के लिए कई एफएफटी के परिणाम को औसत करना चाहते हैं, या क्योंकि आप स्पेक्ट्रा की तुलना करना चाहते हैं संकेतों को बर्बाद किए बिना अलग-अलग अवधि होती है। केवल एनएफएफटी < के मूल्य के साथ एफएफटी कमांड का उपयोग करके आपके सिग्नल में तत्वों की संख्या के परिणामस्वरूप सिग्नल के अंतिम एनएफएफटी बिंदुओं से गणना की गई एफएफटी होगी। यह थोड़ा अपमानजनक है।

निम्नलिखित उदाहरण उपयोगी एप्लिकेशन के लिए अधिक प्रासंगिक है। इससे पता चलता है कि कैसे आप ब्लॉक में एक संकेत विभाजित और फिर प्रत्येक ब्लॉक की प्रक्रिया और परिणाम औसतन:

%//These are the user parameters 
durT = 1; 
fs = 22050; 
NFFT = 2048; 
sigFreq = 300; 

%//Calculate time axis 
dt = 1/fs; 
tAxis = dt:dt:(durT-dt); 

%//Calculate frequency axis 
df = fs/NFFT; 
fAxis = 0:df:(fs-df); 

%//Calculate time domain signal 
x = cos( 2*pi*sigFreq*tAxis ); 

%//Buffer it and window 
win = hamming(NFFT);%//chose window type based on your application 
x = buffer(x, NFFT, NFFT/2); %// 50% overlap between frames in this instance 
x = x(:, 2:end-1); %//optional step to remove zero padded frames 
x = ( x' * diag(win) )'; %//efficiently window each frame using matrix algebra 

%// Calculate mean FFT 
F = abs( fft(x, NFFT)/sum(win) ); 
F = mean(F,2); 

subplot(2,1,1); 
plot( fAxis, 2*F ) 
xlim([0 2*sigFreq]) 
title('single sided spectrum') 

subplot(2,1,2); 
plot( fAxis-fs/2, fftshift(F) ) 
xlim([-2*sigFreq 2*sigFreq]) 
title('whole fft-shifted spectrum') 

मैं ऊपर के उदाहरण में एक आलोचनात्मक खिड़की का उपयोग करें। आपके द्वारा चुनी गई विंडो को http://en.wikipedia.org/wiki/Window_function

आपके द्वारा चुने गए ओवरलैप राशि को आपके द्वारा उपयोग की जाने वाली विंडो के प्रकार पर निर्भर करेगा। उपर्युक्त उदाहरण में, हैमिंग विंडो प्रत्येक बफर में नमूनों को प्रत्येक फ्रेम के केंद्र से शून्य की तरफ दूर करती है। इनपुट सिग्नल में सभी जानकारी का उपयोग करने के लिए, कुछ ओवरलैप का उपयोग करना महत्वपूर्ण है। हालांकि, अगर आप केवल एक सादे आयताकार खिड़की का उपयोग करते हैं, तो ओवरलैप व्यर्थ हो जाता है क्योंकि सभी नमूनों को समान रूप से भारित किया जाता है। जितना अधिक ओवरलैप आप उपयोग करते हैं, उतना अधिक प्रसंस्करण आवश्यक स्पेक्ट्रम की गणना करने के लिए आवश्यक है।

आशा है कि इससे आपकी समझ में मदद मिलेगी।

Answer 2

आपका परिणाम बिल्कुल सही है। आपकी आवृत्ति अक्ष गणना भी सही है। समस्या वाई धुरी पैमाने पर निहित है। जब आप फ़ंक्शन xlims का उपयोग करते हैं, तो मैटलैब स्वचालित रूप से वाई स्केल को पुन: गणना करता है ताकि आप "सार्थक" डेटा देख सकें। जब कोसाइन शिखर आपके द्वारा चुने गए सीमा से बाहर झूठ बोलते हैं (जब एफ> 500 हर्ट्ज), तो दिखाने के लिए कोई चोट नहीं होती है, इसलिए स्केल की गणना कुछ घुटनों के छोटे शोर के आधार पर की जाती है (यहां मेरे कंप्यूटर पर, मैटलैब 2011a के साथ, वाई स्केल 10 था -16)।

सीमा बदलना वास्तव में सही दृष्टिकोण है, क्योंकि यदि आप इसे नहीं बदलते हैं तो आप आवृत्ति स्पेक्ट्रम पर चोटियों को नहीं देख सकते हैं।

हालांकि, मैंने देखा एक बात। क्या आपके पास ट्रांसफॉर्म के वास्तविक हिस्से को साजिश करने का कोई कारण है? आमतौर पर, यह abs(F) है जो प्लॉट हो जाता है, न कि असली हिस्सा।

संपादित करें: वास्तव में, आप आवृत्ति अक्ष केवल सही हैं क्योंकि डीएफ, इस मामले में, 1 है। फैक्सिस लाइन सही है, लेकिन डीएफ गणना नहीं है।

एफएफटी एन अंकों की गणना करता है -एफएस/2 से एफएस/2 तक। तो एफएस की एक श्रृंखला से एन अंक एफएस/एन के एक डीएफ पैदा करता है। एन/समय = एफएस => समय = एन/एफएस के रूप में। डीएफ की अभिव्यक्ति पर आपने इसका उपयोग किया है: your_df = Fs/N * (N/Fs) = (Fs/N)^2। एफएस/एन = 1 के रूप में, अंतिम परिणाम सही था: पी