Notice: This put up is an excerpt from the forthcoming e book, Deep Studying and Scientific Computing with R torch. The chapter in query is on the Discrete Fourier Rework (DFT), and is situated partly three. Half three is devoted to scientific computation past deep studying.
There are two chapters on the Fourier Rework. The primary strives to, in as “verbal” and lucid a method as was potential to me, forged a lightweight on what’s behind the magic; it additionally exhibits how, surprisingly, you possibly can code the DFT in merely half a dozen traces. The second focuses on quick implementation (the Quick Fourier Rework, or FFT), once more with each conceptual/explanatory in addition to sensible, codeityourself elements.
Collectively, these cowl way more materials than might sensibly match right into a weblog put up; due to this fact, please contemplate what follows extra as a “teaser” than a completely fledged article.
Within the sciences, the Fourier Rework is nearly in every single place. Said very usually, it converts knowledge from one illustration to a different, with none lack of info (if executed accurately, that’s.) If you happen to use torch
, it’s only a perform name away: torch_fft_fft()
goes a method, torch_fft_ifft()
the opposite. For the person, that’s handy – you “simply” must know tips on how to interpret the outcomes. Right here, I wish to assist with that. We begin with an instance perform name, enjoying round with its output, after which, attempt to get a grip on what’s going on behind the scenes.
Understanding the output of torch_fft_fft()
As we care about precise understanding, we begin from the only potential instance sign, a pure cosine that performs one revolution over the entire sampling interval.
Start line: A cosine of frequency 1
The way in which we set issues up, there will probably be sixtyfour samples; the sampling interval thus equals N = 64
. The content material of frequency()
, the beneath helper perform used to assemble the sign, displays how we characterize the cosine. Specifically:
[
f(x) = cos(frac{2 pi}{N} k x)
]
Right here (x) values progress over time (or area), and (okay) is the frequency index. A cosine is periodic with interval (2 pi); so if we would like it to first return to its beginning state after sixtyfour samples, and (x) runs between zero and sixtythree, we’ll need (okay) to be equal to (1). Like that, we’ll attain the preliminary state once more at place (x = frac{2 pi}{64} * 1 * 64).
Let’s rapidly affirm this did what it was purported to:
df < knowledge.body(x = sample_positions, y = as.numeric(x))
ggplot(df, aes(x = x, y = y)) +
geom_line() +
xlab("time") +
ylab("amplitude") +
theme_minimal()
Now that we have now the enter sign, torch_fft_fft()
computes for us the Fourier coefficients, that’s, the significance of the assorted frequencies current within the sign. The variety of frequencies thoughtabout will equal the variety of sampling factors: So (X) will probably be of size sixtyfour as properly.
(In our instance, you’ll discover that the second half of coefficients will equal the primary in magnitude. That is the case for each realvalued sign. In such circumstances, you might name torch_fft_rfft()
as a substitute, which yields “nicer” (within the sense of shorter) vectors to work with. Right here although, I wish to clarify the overall case, since that’s what you’ll discover executed in most expositions on the subject.)
Even with the sign being actual, the Fourier coefficients are complicated numbers. There are 4 methods to examine them. The primary is to extract the true half:
[1] 0 32 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
[29] 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
[57] 0 0 0 0 0 0 0 32
Solely a single coefficient is nonzero, the one at place 1. (We begin counting from zero, and should discard the second half, as defined above.)
Now trying on the imaginary half, we discover it’s zero all through:
[1] 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
[29] 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
[57] 0 0 0 0 0 0 0 0
At this level we all know that there’s only a single frequency current within the sign, specifically, that at (okay = 1). This matches (and it higher needed to) the best way we constructed the sign: specifically, as carrying out a single revolution over the entire sampling interval.
Since, in idea, each coefficient might have nonzero actual and imaginary elements, usually what you’d report is the magnitude (the sq. root of the sum of squared actual and imaginary elements):
[1] 0 32 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
[29] 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
[57] 0 0 0 0 0 0 0 32
Unsurprisingly, these values precisely mirror the respective actual elements.
Lastly, there’s the section, indicating a potential shift of the sign (a pure cosine is unshifted). In torch
, we have now torch_angle()
complementing torch_abs()
, however we have to have in mind roundoff error right here. We all know that in every however a single case, the true and imaginary elements are each precisely zero; however on account of finite precision in how numbers are introduced in a pc, the precise values will usually not be zero. As an alternative, they’ll be very small. If we take one among these “faux nonzeroes” and divide it by one other, as occurs within the angle calculation, massive values may end up. To stop this from taking place, our customized implementation rounds each inputs earlier than triggering the division.
section < perform(Ft, threshold = 1e5) {
torch_atan2(
torch_abs(torch_round(Ft$imag * threshold)),
torch_abs(torch_round(Ft$actual * threshold))
)
}
as.numeric(section(Ft)) %>% spherical(5)
[1] 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
[29] 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
[57] 0 0 0 0 0 0 0 0
As anticipated, there is no such thing as a section shift within the sign.
Let’s visualize what we discovered.
create_plot < perform(x, y, amount) {
df < knowledge.body(
x_ = x,
y_ = as.numeric(y) %>% spherical(5)
)
ggplot(df, aes(x = x_, y = y_)) +
geom_col() +
xlab("frequency") +
ylab(amount) +
theme_minimal()
}
p_real < create_plot(
sample_positions,
real_part,
"actual half"
)
p_imag < create_plot(
sample_positions,
imag_part,
"imaginary half"
)
p_magnitude < create_plot(
sample_positions,
magnitude,
"magnitude"
)
p_phase < create_plot(
sample_positions,
section(Ft),
"section"
)
p_real + p_imag + p_magnitude + p_phase
It’s truthful to say that we have now no purpose to doubt what torch_fft_fft()
has executed. However with a pure sinusoid like this, we are able to perceive precisely what’s occurring by computing the DFT ourselves, by hand. Doing this now will considerably assist us later, after we’re writing the code.
Reconstructing the magic
One caveat about this part. With a subject as wealthy because the Fourier Rework, and an viewers who I think about to fluctuate extensively on a dimension of math and sciences schooling, my probabilities to satisfy your expectations, expensive reader, should be very near zero. Nonetheless, I wish to take the danger. If you happen to’re an skilled on this stuff, you’ll anyway be simply scanning the textual content, looking for items of torch
code. If you happen to’re reasonably acquainted with the DFT, you should still like being reminded of its internal workings. And – most significantly – if you happen to’re fairly new, and even utterly new, to this matter, you’ll hopefully take away (not less than) one factor: that what looks as if one of many biggest wonders of the universe (assuming there’s a actuality by some means similar to what goes on in our minds) might be a marvel, however neither “magic” nor a factor reserved to the initiated.
In a nutshell, the Fourier Rework is a foundation transformation. Within the case of the DFT – the Discrete Fourier Rework, the place time and frequency representations each are finite vectors, not capabilities – the brand new foundation appears to be like like this:
[
begin{aligned}
&mathbf{w}^{0n}_N = e^{ifrac{2 pi}{N}* 0 * n} = 1
&mathbf{w}^{1n}_N = e^{ifrac{2 pi}{N}* 1 * n} = e^{ifrac{2 pi}{N} n}
&mathbf{w}^{2n}_N = e^{ifrac{2 pi}{N}* 2 * n} = e^{ifrac{2 pi}{N}2n}& …
&mathbf{w}^{(N1)n}_N = e^{ifrac{2 pi}{N}* (N1) * n} = e^{ifrac{2 pi}{N}(N1)n}
end{aligned}
]
Right here (N), as earlier than, is the variety of samples (64, in our case); thus, there are (N) foundation vectors. With (okay) operating by the idea vectors, they are often written:
[
mathbf{w}^{kn}_N = e^{ifrac{2 pi}{N}k n}
] {#eqdft1}
Like (okay), (n) runs from (0) to (N1). To grasp what these foundation vectors are doing, it’s useful to quickly swap to a shorter sampling interval, (N = 4), say. If we achieve this, we have now 4 foundation vectors: (mathbf{w}^{0n}_N), (mathbf{w}^{1n}_N), (mathbf{w}^{2n}_N), and (mathbf{w}^{3n}_N). The primary one appears to be like like this:
[
mathbf{w}^{0n}_N
=
begin{bmatrix}
e^{ifrac{2 pi}{4}* 0 * 0}
e^{ifrac{2 pi}{4}* 0 * 1}
e^{ifrac{2 pi}{4}* 0 * 2}
e^{ifrac{2 pi}{4}* 0 * 3}
end{bmatrix}
=
begin{bmatrix}
1
1
1
1
end{bmatrix}
]
The second, like so:
[
mathbf{w}^{1n}_N
=
begin{bmatrix}
e^{ifrac{2 pi}{4}* 1 * 0}
e^{ifrac{2 pi}{4}* 1 * 1}
e^{ifrac{2 pi}{4}* 1 * 2}
e^{ifrac{2 pi}{4}* 1 * 3}
end{bmatrix}
=
begin{bmatrix}
1
e^{ifrac{pi}{2}}
e^{i pi}
e^{ifrac{3 pi}{4}}
end{bmatrix}
=
begin{bmatrix}
1
i
1
i
end{bmatrix}
]
That is the third:
[
mathbf{w}^{2n}_N
=
begin{bmatrix}
e^{ifrac{2 pi}{4}* 2 * 0}
e^{ifrac{2 pi}{4}* 2 * 1}
e^{ifrac{2 pi}{4}* 2 * 2}
e^{ifrac{2 pi}{4}* 2 * 3}
end{bmatrix}
=
begin{bmatrix}
1
e^{ipi}
e^{i 2 pi}
e^{ifrac{3 pi}{2}}
end{bmatrix}
=
begin{bmatrix}
1
1
1
1
end{bmatrix}
]
And eventually, the fourth:
[
mathbf{w}^{3n}_N
=
begin{bmatrix}
e^{ifrac{2 pi}{4}* 3 * 0}
e^{ifrac{2 pi}{4}* 3 * 1}
e^{ifrac{2 pi}{4}* 3 * 2}
e^{ifrac{2 pi}{4}* 3 * 3}
end{bmatrix}
=
begin{bmatrix}
1
e^{ifrac{3 pi}{2}}
e^{i 3 pi}
e^{ifrac{9 pi}{2}}
end{bmatrix}
=
begin{bmatrix}
1
i
1
i
end{bmatrix}
]
We will characterize these 4 foundation vectors when it comes to their “pace”: how briskly they transfer across the unit circle. To do that, we merely have a look at the rightmost column vectors, the place the ultimate calculation outcomes seem. The values in that column correspond to positions pointed to by the revolving foundation vector at totally different cutoff dates. Because of this a single “replace of place”, we are able to see how briskly the vector is shifting in a single time step.
Wanting first at (mathbf{w}^{0n}_N), we see that it doesn’t transfer in any respect. (mathbf{w}^{1n}_N) goes from (1) to (i) to (1) to (i); yet another step, and it might be again the place it began. That’s one revolution in 4 steps, or a step measurement of (frac{pi}{2}). Then (mathbf{w}^{2n}_N) goes at double that tempo, shifting a distance of (pi) alongside the circle. That method, it finally ends up finishing two revolutions general. Lastly, (mathbf{w}^{3n}_N) achieves three full loops, for a step measurement of (frac{3 pi}{2}).
The factor that makes these foundation vectors so helpful is that they’re mutually orthogonal. That’s, their dot product is zero:
[
langle mathbf{w}^{kn}_N, mathbf{w}^{ln}_N rangle = sum_{n=0}^{N1} ({e^{ifrac{2 pi}{N}k n}})^* e^{ifrac{2 pi}{N}l n} = sum_{n=0}^{N1} ({e^{ifrac{2 pi}{N}k n}})e^{ifrac{2 pi}{N}l n} = 0
] {#eqdft2}
Let’s take, for instance, (mathbf{w}^{2n}_N) and (mathbf{w}^{3n}_N). Certainly, their dot product evaluates to zero.
[
begin{bmatrix}
1 & 1 & 1 & 1
end{bmatrix}
begin{bmatrix}
1
i
1
i
end{bmatrix}
=
1 + i + (1) + (i) = 0
]
Now, we’re about to see how the orthogonality of the Fourier foundation considerably simplifies the calculation of the DFT. Did you discover the similarity between these foundation vectors and the best way we wrote the instance sign? Right here it’s once more:
[
f(x) = cos(frac{2 pi}{N} k x)
]
If we handle to characterize this perform when it comes to the idea vectors (mathbf{w}^{kn}_N = e^{ifrac{2 pi}{N}okay n}), the internal product between the perform and every foundation vector will probably be both zero (the “default”) or a a number of of 1 (in case the perform has a part matching the idea vector in query). Fortunately, sines and cosines can simply be transformed into complicated exponentials. In our instance, that is how that goes:
[
begin{aligned}
mathbf{x}_n &= cos(frac{2 pi}{64} n)
&= frac{1}{2} (e^{ifrac{2 pi}{64} n} + e^{ifrac{2 pi}{64} n})
&= frac{1}{2} (e^{ifrac{2 pi}{64} n} + e^{ifrac{2 pi}{64} 63n})
&= frac{1}{2} (mathbf{w}^{1n}_N + mathbf{w}^{63n}_N)
end{aligned}
]
Right here step one instantly outcomes from Euler’s formulation, and the second displays the truth that the Fourier coefficients are periodic, with frequency 1 being the identical as 63, 2 equaling 62, and so forth.
Now, the (okay)th Fourier coefficient is obtained by projecting the sign onto foundation vector (okay).
As a result of orthogonality of the idea vectors, solely two coefficients is not going to be zero: these for (mathbf{w}^{1n}_N) and (mathbf{w}^{63n}_N). They’re obtained by computing the internal product between the perform and the idea vector in query, that’s, by summing over (n). For every (n) ranging between (0) and (N1), we have now a contribution of (frac{1}{2}), leaving us with a closing sum of (32) for each coefficients. For instance, for (mathbf{w}^{1n}_N):
[
begin{aligned}
X_1 &= langle mathbf{w}^{1n}_N, mathbf{x}_n rangle
&= langle mathbf{w}^{1n}_N, frac{1}{2} (mathbf{w}^{1n}_N + mathbf{w}^{63n}_N) rangle
&= frac{1}{2} * 64
&= 32
end{aligned}
]
And analogously for (X_{63}).
Now, trying again at what torch_fft_fft()
gave us, we see we have been in a position to arrive on the identical end result. And we’ve discovered one thing alongside the best way.
So long as we stick with alerts composed of a number of foundation vectors, we are able to compute the DFT on this method. On the finish of the chapter, we’ll develop code that can work for all alerts, however first, let’s see if we are able to dive even deeper into the workings of the DFT. Three issues we’ll wish to discover:

What would occur if frequencies modified – say, a melody have been sung at the next pitch?

What about amplitude modifications – say, the music have been performed twice as loud?

What about section – e.g., there have been an offset earlier than the piece began?
In all circumstances, we’ll name torch_fft_fft()
solely as soon as we’ve decided the end result ourselves.
And eventually, we’ll see how complicated sinusoids, made up of various parts, can nonetheless be analyzed on this method, offered they are often expressed when it comes to the frequencies that make up the idea.
Various frequency
Assume we quadrupled the frequency, giving us a sign that seemed like this:
[
mathbf{x}_n = cos(frac{2 pi}{N}*4*n)
]
Following the identical logic as above, we are able to categorical it like so:
[
mathbf{x}_n = frac{1}{2} (mathbf{w}^{4n}_N + mathbf{w}^{60n}_N)
]
We already see that nonzero coefficients will probably be obtained just for frequency indices (4) and (60). Choosing the previous, we acquire
[
begin{aligned}
X_4 &= langle mathbf{w}^{4n}_N, mathbf{x}_n rangle
&= langle mathbf{w}^{4n}_N, frac{1}{2} (mathbf{w}^{4n}_N + mathbf{w}^{60n}_N) rangle
&= 32
end{aligned}
]
For the latter, we’d arrive on the identical end result.
Now, let’s be certain our evaluation is appropriate. The next code snippet accommodates nothing new; it generates the sign, calculates the DFT, and plots them each.
x < torch_cos(frequency(4, N) * sample_positions)
plot_ft < perform(x) p_imag) /
(p_magnitude
plot_ft(x)
This does certainly affirm our calculations.
A particular case arises when sign frequency rises to the very best one “allowed”, within the sense of being detectable with out aliasing. That would be the case at one half of the variety of sampling factors. Then, the sign will seem like so:
[
mathbf{x}_n = frac{1}{2} (mathbf{w}^{32n}_N + mathbf{w}^{32n}_N)
]
Consequently, we find yourself with a single coefficient, similar to a frequency of 32 revolutions per pattern interval, of double the magnitude (64, thus). Listed here are the sign and its DFT:
x < torch_cos(frequency(32, N) * sample_positions)
plot_ft(x)
Various amplitude
Now, let’s take into consideration what occurs after we fluctuate amplitude. For instance, say the sign will get twice as loud. Now, there will probably be a multiplier of two that may be taken exterior the internal product. In consequence, the one factor that modifications is the magnitude of the coefficients.
Let’s confirm this. The modification is predicated on the instance we had earlier than the final one, with 4 revolutions over the sampling interval:
x < 2 * torch_cos(frequency(4, N) * sample_positions)
plot_ft(x)
To this point, we have now not as soon as seen a coefficient with nonzero imaginary half. To alter this, we add in section.
Including section
Altering the section of a sign means shifting it in time. Our instance sign is a cosine, a perform whose worth is 1 at (t=0). (That additionally was the – arbitrarily chosen – place to begin of the sign.)
Now assume we shift the sign ahead by (frac{pi}{2}). Then the height we have been seeing at zero strikes over to (frac{pi}{2}); and if we nonetheless begin “recording” at zero, we should discover a worth of zero there. An equation describing that is the next. For comfort, we assume a sampling interval of (2 pi) and (okay=1), in order that the instance is a straightforward cosine:
[
f(x) = cos(x – phi)
]
The minus signal might look unintuitive at first. But it surely does make sense: We now wish to acquire a worth of 1 at (x=frac{pi}{2}), so (x – phi) ought to consider to zero. (Or to any a number of of (pi).) Summing up, a delay in time will seem as a unfavourable section shift.
Now, we’re going to calculate the DFT for a shifted model of our instance sign. However if you happen to like, take a peek on the phaseshifted model of the timedomain image now already. You’ll see {that a} cosine, delayed by (frac{pi}{2}), is nothing else than a sine beginning at 0.
To compute the DFT, we comply with our familiarbynow technique. The sign now appears to be like like this:
[
mathbf{x}_n = cos(frac{2 pi}{N}*4*x – frac{pi}{2})
]
First, we categorical it when it comes to foundation vectors:
[
begin{aligned}
mathbf{x}_n &= cos(frac{2 pi}{64} 4 n – frac{pi}{2})
&= frac{1}{2} (e^{ifrac{2 pi}{64} 4n – frac{pi}{2}} + e^{ifrac{2 pi}{64} 60n – frac{pi}{2}})
&= frac{1}{2} (e^{ifrac{2 pi}{64} 4n} e^{i frac{pi}{2}} + e^{ifrac{2 pi}{64} 60n} e^{ifrac{pi}{2}})
&= frac{1}{2} (e^{i frac{pi}{2}} mathbf{w}^{4n}_N + e^{i frac{pi}{2}} mathbf{w}^{60n}_N)
end{aligned}
]
Once more, we have now nonzero coefficients just for frequencies (4) and (60). However they’re complicated now, and each coefficients are not equivalent. As an alternative, one is the complicated conjugate of the opposite. First, (X_4):
[
begin{aligned}
X_4 &= langle mathbf{w}^{4n}_N, mathbf{x}_n rangle
&=langle mathbf{w}^{4n}_N, frac{1}{2} (e^{i frac{pi}{2}} mathbf{w}^{4n}_N + e^{i frac{pi}{2}} mathbf{w}^{60n}_N) rangle
&= 32 *e^{i frac{pi}{2}}
&= 32i
end{aligned}
]
And right here, (X_{60}):
[
begin{aligned}
X_{60} &= langle mathbf{w}^{60n}_N, mathbf{x}_N rangle
&= 32 *e^{i frac{pi}{2}}
&= 32i
end{aligned}
]
As typical, we examine our calculation utilizing torch_fft_fft()
.
x < torch_cos(frequency(4, N) * sample_positions  pi / 2)
plot_ft(x)
For a pure sine wave, the nonzero Fourier coefficients are imaginary. The section shift within the coefficients, reported as (frac{pi}{2}), displays the time delay we utilized to the sign.
Lastly – earlier than we write some code – let’s put all of it collectively, and have a look at a wave that has greater than a single sinusoidal part.
Superposition of sinusoids
The sign we assemble should still be expressed when it comes to the idea vectors, however it’s not a pure sinusoid. As an alternative, it’s a linear mixture of such:
[
begin{aligned}
mathbf{x}_n &= 3 sin(frac{2 pi}{64} 4n) + 6 cos(frac{2 pi}{64} 2n) +2cos(frac{2 pi}{64} 8n)
end{aligned}
]
I gained’t undergo the calculation intimately, however it’s no totally different from the earlier ones. You compute the DFT for every of the three parts, and assemble the outcomes. With none calculation, nevertheless, there’s fairly just a few issues we are able to say:
 For the reason that sign consists of two pure cosines and one pure sine, there will probably be 4 coefficients with nonzero actual elements, and two with nonzero imaginary elements. The latter will probably be complicated conjugates of one another.
 From the best way the sign is written, it’s straightforward to find the respective frequencies, as properly: The allreal coefficients will correspond to frequency indices 2, 8, 56, and 62; the allimaginary ones to indices 4 and 60.
 Lastly, amplitudes will end result from multiplying with (frac{64}{2}) the scaling components obtained for the person sinusoids.
Let’s examine:
Now, how will we calculate the DFT for much less handy alerts?
Coding the DFT
Happily, we already know what needs to be executed. We wish to mission the sign onto every of the idea vectors. In different phrases, we’ll be computing a bunch of internal merchandise. Logicwise, nothing modifications: The one distinction is that usually, it is not going to be potential to characterize the sign when it comes to just some foundation vectors, like we did earlier than. Thus, all projections will really need to be calculated. However isn’t automation of tedious duties one factor we have now computer systems for?
Let’s begin by stating enter, output, and central logic of the algorithm to be applied. As all through this chapter, we keep in a single dimension. The enter, thus, is a onedimensional tensor, encoding a sign. The output is a onedimensional vector of Fourier coefficients, of the identical size because the enter, every holding details about a frequency. The central thought is: To acquire a coefficient, mission the sign onto the corresponding foundation vector.
To implement that concept, we have to create the idea vectors, and for every one, compute its internal product with the sign. This may be executed in a loop. Surprisingly little code is required to perform the aim:
dft < perform(x) {
n_samples < size(x)
n < torch_arange(0, n_samples  1)$unsqueeze(1)
Ft < torch_complex(
torch_zeros(n_samples), torch_zeros(n_samples)
)
for (okay in 0:(n_samples  1)) {
w_k < torch_exp(1i * 2 * pi / n_samples * okay * n)
dot < torch_matmul(w_k, x$to(dtype = torch_cfloat()))
Ft[k + 1] < dot
}
Ft
}
To check the implementation, we are able to take the final sign we analysed, and evaluate with the output of torch_fft_fft()
.
[1] 0 0 192 0 0 0 0 0 64 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
[29] 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
[57] 64 0 0 0 0 0 192 0
[1] 0 0 0 0 96 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
[29] 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
[57] 0 0 0 0 96 0 0 0
Reassuringly – if you happen to look again – the outcomes are the identical.
Above, did I say “little code”? Actually, a loop isn’t even wanted. As an alternative of working with the idea vectors onebyone, we are able to stack them in a matrix. Then every row will maintain the conjugate of a foundation vector, and there will probably be (N) of them. The columns correspond to positions (0) to (N1); there will probably be (N) of them as properly. For instance, that is how the matrix would search for (N=4):
[
mathbf{W}_4
=
begin{bmatrix}
e^{ifrac{2 pi}{4}* 0 * 0} & e^{ifrac{2 pi}{4}* 0 * 1} & e^{ifrac{2 pi}{4}* 0 * 2} & e^{ifrac{2 pi}{4}* 0 * 3}
e^{ifrac{2 pi}{4}* 1 * 0} & e^{ifrac{2 pi}{4}* 1 * 1} & e^{ifrac{2 pi}{4}* 1 * 2} & e^{ifrac{2 pi}{4}* 1 * 3}
e^{ifrac{2 pi}{4}* 2 * 0} & e^{ifrac{2 pi}{4}* 2 * 1} & e^{ifrac{2 pi}{4}* 2 * 2} & e^{ifrac{2 pi}{4}* 2 * 3}
e^{ifrac{2 pi}{4}* 3 * 0} & e^{ifrac{2 pi}{4}* 3 * 1} & e^{ifrac{2 pi}{4}* 3 * 2} & e^{ifrac{2 pi}{4}* 3 * 3}
end{bmatrix}
] {#eqdft3}
Or, evaluating the expressions:
[
mathbf{W}_4
=
begin{bmatrix}
1 & 1 & 1 & 1
1 & i & 1 & i
1 & 1 & 1 & 1
1 & i & 1 & i
end{bmatrix}
]
With that modification, the code appears to be like much more elegant:
dft_vec < perform(x) {
n_samples < size(x)
n < torch_arange(0, n_samples  1)$unsqueeze(1)
okay < torch_arange(0, n_samples  1)$unsqueeze(2)
mat_k_m < torch_exp(1i * 2 * pi / n_samples * okay * n)
torch_matmul(mat_k_m, x$to(dtype = torch_cfloat()))
}
As you possibly can simply confirm, the end result is identical.
Thanks for studying!
Picture by Trac Vu on Unsplash