필터 지우기
필터 지우기

FFT of NMR signal: frequency cut off

조회 수: 15 (최근 30일)
Nicolai
Nicolai 2014년 6월 28일
답변: Aaron Globisch 2017년 5월 15일
Hello,
I've got a problem with the fourier transform of a time domain nmr signal. When I perform the fourier transform I receive a weired spectrum. The first Peak is cutted in half, as you can see in the attached picture. Does somebody know, why lower frequencies are not preserved? Is there a mistake in my code or does somebody have an idea? The raw data is ok and contain all informations as I know when I process the data with standard nmr software.
Many thanks, N.
clear all;
close all;
clc
tdomain = load('nmr_ethanol.mat');
tdomain = tdomain.signal;
AqT = tdomain(1,end);
SampRate = tdomain(2,1)-tdomain(1,1);
maxfreq = (length(tdomain(:,1))-1)./(length(tdomain(:,1))*SampRate);
lowfrq = -2295.07256526605;
FreqVek = linspace(0,maxfreq,length(tdomain(:,1)))+(lowfrq);
FreqVek2 = linspace(0,maxfreq,2*length(tdomain(:,1)))+(lowfrq);
npoints = 2.*length(tdomain(:,1));
% fourier transformation
fdomain = fft(tdomain(:,2),npoints);
% cut off symmetrical part and frequency shift
fdomain_sym = fdomain(1:npoints./2);
fdomain_shift = fftshift(fdomain_sym/npoints);
% Plot
fig1 = figure(1);
subplot(2,2,1)
plot(tdomain(:,1),tdomain(:,2));
title('time domain')
xlabel('time, t / s');
xlim([0 0.75])
subplot(2,2,2)
plot(FreqVek2,real(fdomain),'-');
grid on;
title('frequency domain')
xlabel('frequency, f / s^{-1}');
xlim([-2500 3300]);
subplot(2,2,3)
plot(FreqVek,real(fdomain_sym),'-');
grid on;
title('frequency domain, symmetrical half')
xlabel('frequency, f / s^{-1}');
xlim([-2500 -1500]);
subplot(2,2,4)
plot(FreqVek,real(fdomain_shift),'-');
grid on;
title('frequency domain, frequency shift')
xlabel('frequency, f / s^{-1}');
xlim([200 1000]);

답변 (4개)

Nicolai
Nicolai 2014년 6월 29일
Thanks a lot for your answers. Unfortunatly the change in the plotting settings do not solve the problem. In nmr spectroscopy we always look at the real part of the spectrum and do not use the magnitude spectrum. But thanks for the suggestions!
Could there be a reason why I cannot restore all frequencies?
The spectrum supposed to be look like this:
Thanks a lot for all your ideas! N.
  댓글 수: 1
Aaron Globisch
Aaron Globisch 2017년 5월 10일
Hi Nicolai,
did you find a solution for your problem?

댓글을 달려면 로그인하십시오.


Star Strider
Star Strider 2014년 6월 28일
Not certain about this because I don’t have your data, but see if changing your plot statements in the last two supblot calls from real to abs:
plot(FreqVek,abs(fdomain_sym),'-');
and
plot(FreqVek,abs(fdomain_shift),'-');
solves your problem.

Daniel kiracofe
Daniel kiracofe 2014년 6월 28일
Hi Nicolai. I'm not sure exactly what you are expecting to see, but I suspect you just didn't multiply by 2 when converting from 2 sided to 1 sided spectrum. Generally, I never use matlab's fft() directly. I always use a wrapper script. Below is a basic wrapper script. You can find more discussion of it at my website <http://www.mechanicalvibration.com/Making_matlab_s_fft_functio.html>
% function [frq, amp, phase] = simpleFFT( signal, ScanRate)
% Purpose: perform an FFT of a real-valued input signal, and generate the single-sided
% output, in amplitude and phase, scaled to the same units as the input.
%inputs:
% signal: the signal to transform
% ScanRate: the sampling frequency (in Hertz)
% outputs:
% frq: a vector of frequency points (in Hertz)
% amp: a vector of amplitudes (same units as the input signal)
% phase: a vector of phases (in radians)
function [frq, amp, phase] = simpleFFT( signal, ScanRate)
n = length(signal);
z = fft(signal, n); %do the actual work
%generate the vector of frequencies
halfn = n / 2;
deltaf = 1 / ( n / ScanRate);
frq = (0:(halfn-1)) * deltaf;
% convert from 2 sided spectrum to 1 sided
%(assuming that the input is a real signal)
amp(1) = abs(z(1)) ./ (n);
amp(2:halfn) = abs(z(2:halfn)) ./ (n / 2);
phase = angle(z(1:halfn));

Aaron Globisch
Aaron Globisch 2017년 5월 15일
Hello,
i figured out the correct way:
specturm = abs(fftshift(fft(data, s_points)/s_points));
plot(FreqVek2, spectrum);
where s_points is the amount of sampled points and FreqVek2 the frequency vector just like u declared it.

카테고리

Help CenterFile Exchange에서 Parametric Spectral Estimation에 대해 자세히 알아보기

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by