Tag Info

Hot answers tagged

23

I think there are at least three elements to consider here: FourierTransform and Fourier, by default, output results in different forms Plotting Sin[x] UnitStep[x] is not the same as Sin[x] and behaves differently when used in conjunction with Fourier and FourierTransform Plot does not handle DiracDelta elegantly The signal processing form of the ...


18

This is how I would do this. Define frequencies and sampling rate precisely. Then use Periodogram because it takes SampleRate as an option and rescales frequency axis automatically. Read up Docs on Periodogram - see examples there. data = Table[{t, Sin[2 Pi 697 t] + Sin[2 Pi 1209 t]}, {t, 0., 0.1, 1/8000.}]; ListLinePlot[data, AspectRatio -> 1/4, ...


17

Based on the MATLAB documentation, it would appear that this is accomplished by simple zero-filling. As such, you can obtain the same result in Mathematica using Fourier[PadRight[list, n, 0.], FourierParameters -> {1, -1}] where list is your signal and n is the desired length. For a multidimensional FFT, replace n with {n1, n2, ...}, where n1, n2, ...


16

Here's a possible starting point for a solution. It splits the sample list into chunks and measures the Norm of the sample Differences in each chunk, and then does the FFT on that data. bpmplot[snd_, bpmmax_: 300] := Module[{samples, minfreq, signal, fft}, samples = snd[[1, 1, 1]]; minfreq = snd[[1, 2]]/Length[samples]; signal = (Norm[Differences[#]]) ...


15

There is the function NFourierTransform[] (as well as NInverseFourierTransform[]) implemented in the package FourierSeries`. The function, as with the related kernel functions, takes a FourierParameters option so you can adjust computations to your preferred normalization as needed. For your specific normalization, you apparently want the setting ...


14

The problem here is that Mathematica doesn't recognize {x, y, z} as some kind of a vector object that should be treated as grouped together; instead, it substitutes in three independent variables, and probably starts integrating them one by one. The result is a very complicated integral. If you do the coordinate transformation yourself, you can reproduce ...


13

I think perhaps you need codes like this: Func[x_] := Sin[x]; tmin = 0; tmax = 10; \[CapitalDelta]t = (tmax - tmin)/100; tgrid = Table[t, {t, tmin, tmax, \[CapitalDelta]t}]; wgrid = RotateRight[(2 \[Pi])/(tmax - tmin)* Range[-((Length@tgrid - 1)/2), (Length@tgrid - 1)/2], ( Length@tgrid - 1)/2]; ListLogLogPlot[{wgrid, (tmax - tmin)/Sqrt[2 ...


12

It's not a bug, it's a feature Exact integration returns 1/Sqrt[2 Pi] Integrate[(1 + Sqrt[x])^2 Exp[I k x], {x, -Infinity, Infinity}, Assumptions -> {k \[Element] Reals}] Integrate::idiv: "Integral of E^(I\k\x)\ (1+[Sqrt]x)^2 does not converge on {-Infinity,Infinity}." However we can multiply by Exp[-b Abs[x]] and then put b -> 0 ...


12

The Fourier transform is defined as: $$ H(f)=\int_\infty^\infty h(t) e^{2\pi i f t}dt\\ h(t)=\int_\infty^\infty H(f) e^{-2\pi i f t}df $$ where $h(t)$ is the signal, and $H(f)$ is it's Fourier transform, if $t$ is meassured in second, then $f$ is measured in Hz. The discrete Fourier transform is defined as: $$ H_{f_j}=\frac{1}{N}\sum_{k}h_{t_k}e^{2\pi i ...


12

Looking at your plotted data you can see about 40 cycles of the dominant frequency, this tells you that the peak will appear somewhere around the 40th element of the DFT. That's in the region where your plot of the DFT is clipped, so it's no wonder you can't see the peak. Looking at the relevant part of the DFT you can see the peak quite clearly: ...


11

In Mathematica there is a designated function for this, UnitBox, PiecewiseExpand[UnitBox[x]] which gives expected result without assumptions: FourierTransform[UnitBox[x/a], x, k, FourierParameters -> {1, -1}] Abs[a] Sinc[(a k)/2] There is actually a set of designated functions: FourierTransform[UnitTriangle[x/a], x, t, FourierParameters ...


11

In the definition of s you're summing from k==0. Since the summand has a term 1/k this gives a divide-by-zero error when calculating the partial sums. The sum should in fact start from k==1 (the zeroth coefficient is taken care of by the constant term in front of the sum). The first few approximations then look like s[n_, x_] := 8/4 + 3/(9 \[Pi]) Sum[(6 ...


11

You could simply remove the vertical waves (e.g. by subtracting the median of each column) and histogram modification. Using @bill s's cropped image: img = Image[ ImageData[ ColorConvert[Import["http://i.stack.imgur.com/EvjuW.png"], "Grayscale"]][[;; , ;; , 1]]]; (* remove alpha channel *) columnMedian = Median[ImageData[img]]; medianRemoved = # ...


10

I also would expect Mathematica to simplify all Fourier transformed derivatives equally, but it may be understandable that the simplifications are harder to see when the derivative is not taken with respect to the innermost Fourier transform variable. To work around this problem, you could change the order of integrations for the Fourier transform to ...


10

Assuming sample rate is 2.8 second per sample, then may be this: SetDirectory[NotebookDirectory[]]; data = Import["testdata.txt", "List"]; ListPlot[data, Joined -> True] py = Fourier[data, FourierParameters -> {1, -1}]; nSamples = Length[data]; nUniquePts = Ceiling[(nSamples + 1)/2]; py = py[[1 ;; nUniquePts]]; py = Abs[py]; py[[1]] = 0; (*zero ...


10

Paying close attention to the documentation for Fourier and FourierTransform one notes that the coefficients of the Sum/Integral terms are different; therefore, to obtain a discrete transform with amplitudes equal to those from the continuous transform, one must multiply the former by Sqrt[2 Pi / n] where n is the length of the dataset: The continuous ...


9

I feel there may be a few issues here. First, you're using FourierDST, the discrete sine transform. I'm not too familiar with this one, but it looks like you shouldn't confuse it with Fourier. Application of FourierDST as follows: ListLinePlot[ FourierDST[Table[Sin[100 t], {t, 0, 10, 0.02}]][[250 ;; 350]], PlotRange -> All] yields: whereas, with ...


9

I had a play with various Compile options and didn't get anywhere (I managed to make it slower though!). However, you can get a nice little speed boost using ParallelTable. Your original on my machine: NFourierTransform[f_Function, {kmin_, kmax_}] := Interpolation@ Table[{k, Chop@NIntegrate[f@x E^(-I k x), {x, -Infinity, Infinity}]}, {k, ...


9

What is happening in your second example (with the single sine wave giving the "beating") is that you have exceeded the Nyquist frequency: what you are seeing is called aliasing. Here's a simple way to explore this (using your DFT function): dt = 0.66125; f[w_, x_] := Sin[w x]; Manipulate[ls = Table[f[w, x], {x, dt, 200 dt, dt}]; Column[{ListPlot[ls, ...


9

update Just to clean things up a bit, we can use the discussion here to make a couple functions that help extract the frequency data from this dataset. I define two functions findPeriod and reconstruct: Clear[findPeriod]; findPeriod[data_, threshold_] := Module[{fs, s1, s = {}, i, a0f, af, pf, pos, fr, frpos, fdata, fdatac, n, per}, n = ...


8

As mentioned by @Rahul, you have not sampled your sine wave often enough and have introduced artifacts due to aliasing. The frequency of Sin[500 x]=Sin[2 Pi f x] is $f=500/(2\pi)$, which is about 80 Hz. At least two samples per cycle are required to avoid aliasing, hence the default $x$ interval of 1 in {x,0,100} must be reduced to less than about ...


8

There is no inconsistency here. If you assume a to be positive, as Daniel mentioned, you get the same answer: FourierTransform[UnitStep[a/2 + x] UnitStep[a/2 - x], x, k, FourierParameters -> {1, -1}, Assumptions -> a > 0] (* (2 Sin[(a k)/2])/k *) Integrate[UnitStep[a/2 + x] UnitStep[a/2 - x] Exp[-I k x], {x, -∞, ∞}, Assumptions -> a ...


8

If you take a look at the documentation, Mathematica's symbolic Fourier transform function, FourierTransform, computes $$\hat f(k) = \frac{1}{\sqrt{2\pi}}\int_{-\infty}^{\infty}f(x)e^{ikx}\mathrm{d}x$$ You can discretize some piece of this integral by limiting $x$ and $k$ to values $x_1 + (r-1)\Delta x$ and $(s-1)\Delta k$ respectively, where $\Delta ...


8

This is borrowed from comp.soft-sys.math.mathematica posts, primarily by Szesi Mukasa. The outside factor I cribbed from a text book. fft[ll_] := Exp[I*Pi*(Range[Length[ll]] + Boole[OddQ[Length[ll]]])]* RotateRight[Fourier[ll, FourierParameters -> {0, 1}], Quotient[Length[ll], 2]] For an example, you could use a Gaussian. f[x_] := Exp[-x^2] ...


8

Something strange is going on here. Here is a computation which illustrates the issue without some of the extraneous aspects. wrong = FourierCoefficient[1/(x^2 + 1), x, 1] The variable wrong now contains what FourierCoefficient thinks is the coefficient of $e^{i x}$ in the fourier series of $1/(x^2+1)$. According to the documentation for ...


8

Finally I found the most promising algorithm proposed in this really good reference Manuel Guizar-Sicairos and Julio C. Gutiérrez-Vega, "Computation of quasi-discrete Hankel transforms of integer order for propagating optical wave fields," J. Opt. Soc. Am. A 21, 53-58 (2004). The authors call the algorithm pth-order quasi-discrete Hankel Transform (pQDHT) ...


8

img = Import["ExampleData/lena.tif"]; Image[img, ImageSize -> 300] data = ImageData[img];(*get data*) {nRow, nCol, nChannel} = Dimensions[data]; d = data[[All, All, 2]]; d = d*(-1)^Table[i + j, {i, nRow}, {j, nCol}]; fw = Fourier[d, FourierParameters -> {1, 1}]; (*adjust for better viewing as needed*) fudgeFactor = 100; abs = fudgeFactor*Log[1 + ...


7

You can use FourierCoefficient to pre-calculate the Fourier coefficients to arbitrary degree and then use the result very effectively. There may be some issues with zero-th degree, therefore I excluded this using Piecewise. Here's the main code block: f[x_] := Piecewise[{ {-x^3 - 2 x, -2 < x < 0}, {-1 + x, 0 <= x <= 2}}, 0 ]; ...


7

Most of your questions can be answered experimentally. You can find out a lot about Mathematica by just interactively playing and timing results. Let's see how it works out in this case: does it automatically optimize for the type of input data (just integers 0 and 1)? d1 = RandomReal[1, 10^7]; Fourier[d1]; // AbsoluteTiming // First (* ==> ...


7

Some deconvolve algorithm can not recover fine details very well, which is of course a shortcoming. However, by taking advantage of this kind of shortcoming, we can separate fine detailed structures from large scale structures effectively. We borrow the same image loader from nikie's answer: img = Image[ ImageData[ColorConvert[ ...



Only top voted, non community-wiki answers of a minimum length are eligible