FFT of NMR signal: frequency cut off

10 次查看（过去 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.
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 Center 和 File Exchange 中查找有关 Spectral Measurements 的更多信息

标签

产品

Community Treasure Hunt

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

Start Hunting!

Translated by