Hacker Newsnew | past | comments | ask | show | jobs | submitlogin

There’s an alternate way to save space: the discrete Hartley transform.

https://news.ycombinator.com/item?id=27386319

https://twitter.com/theshawwn/status/1400383554673065984

  dht(x) = fft( (1 + 1j)*x ).real / area(x)**0.5

  area(x) = np.prod(np.shape(x))
The dht is its own inverse, which is a neat trick. No need for an inverse fft. In other words, dht(dht(x)) == x.

Notice that the output is the real component of the fft. That means the storage required is exactly identical to the input signal. And since the input signal is copied to the real and imaginary components (1 + 1j)*x, the intermediate representation size is also equal to the input signal size. (Just take sin and cos of the input signal.)

The dht has almost all the same properties of the fft. You can blur an image by taking a hamming window and multiplying. You can resize an image by fftshift + pad with zeros + fftshift. Etc.

In practice the dht seems to have all the benefits of fft with very few downsides. I’m amazed it isn’t more popular.



The convolution theorem is much uglier for the hartley transform. I guess this is the main reason. You can see the Fourier transform as a simple variant of the Hartley transform where the linear filtering can be done directly in the frequency domain, without all the fuss of the hartley transform convolutions.

Another advantage of Fourier with respect to Harley is that derivatives become products with polynomials. Thus the Fourier transform transforms linear pde into trivial algebraic equations. For the Hartley transform, you get more complicated functional equations, so that it's not as useful.

> You can blur an image by taking a hamming window and multiplying.

It doesn't seem like you can? (Unless your image has some symmetries).


Here's an image: https://i.imgur.com/hNCwisy.png

Here's fftshift(dht2(img)): https://i.imgur.com/EB6yIb4.png

Here's a hamming window: https://i.imgur.com/NsGKb0X.png

Here's a hamming window raised to power of 20: https://i.imgur.com/ejCviD7.png

Here's the result of multiplying that with the transformed image: https://i.imgur.com/trHdWZS.png

And here's dht(fftshift(that)): https://i.imgur.com/of4bvhv.png

Presto, blurry image. No symmetries required.

We can increase the blur factor by changing the power from 20 to 200: https://i.imgur.com/yqkR3qb.png

And 2000 (megablur): https://i.imgur.com/7h9patv.png

The code I used to blur:

    def blur(img, factor):
        x = dht2(img)
        x = fftshift(x)
        x = x * hamming(img) ** factor
        x = fftshift(x)
        x = dht2(x)
        return x
There doesn't seem to be any limitations -- anything you can do with FFT, you can do here. And since the DHT isn't complex, you don't need to deal with phase/amplitude.

It's kind of like working with wavelet coefficients in a jpg transform, except wavelets are a lot harder to work with compared to fft signals (which is why people use fft instead of wavelets everywhere). But dht is the best of both worlds -- all the power of fft, none of the hassle of complex values.


> There doesn't seem to be any limitations -- anything you can do with FFT, you can do here.

There's certainly no limitations if all you want to do is to convolve your image with symmetric kernels. But you asked for advantages of the Fourier transform, if any. Convolution of two non-symmetric images is one of them ;)

There's also the fact that derivatives have an easier relationship with the DFT than the DHT. If you are interested in discrete linear pde, that's a big advantage.

For signal processing, the main advantage of the DFT is that the amplitudes and phases are cleanly separated. For the applications where you only want the amplitude, it is already there in the DFT, but you have to untangle it from the DHT, because the DHT intermingles the phases and amplitudes.

An example: imagine that you have an image affected by periodic noise, and you want to remove this noise. This periodic noise appears in the DFT as a frequency with very large amplitude (and an arbitrary phase, corresponding to the phase of the periodic noise, that you don't really care about). Thus, it is easy to locate and filter-out this periodic noise using the DFT: just find the frequency with large amplitude and set it to zero.

Now, for the DHT, the frequency that corresponds to the noise period is not necessarily large (in magnitude). It may even be zero, depending on the phase of the noise! To find that frequency, you have to essentially recover the DFT and compute the amplitudes.

Of course the DHT and the DFT are equivalent and there's an easy conversion between both, so you can do exactly the same things with either. The difference is just a matter of convenience. It's certainly useful to know both.

If you are scared of complex numbers for some reason (e.g., if your programming language makes it difficult to work with them) then you may prefer the DHT. Otherwise, the DFT neatly separates the phases and amplitudes and it's slightly easier to interpret.


On this point, for images much of the structure is encoded in the phase. This fact was very surprising to me when I learned it.

See the images in this stack overflow for an example [0].

[0] https://stackoverflow.com/a/40387839


Good morning!

> If you are scared of complex numbers for some reason (e.g., if your programming language makes it difficult to work with them) then you may prefer the DHT.

It’s not about that. The DHT cuts memory bandwidth in half. With FFT it’s a complex number, so it’s doubled.

FFT on a real signal can be cut in half too, but it requires doing weird shuffling of the values. (See submission’s blog post on the topic.) It’s unnecessary complexity unless you really need a full FFT for some very specific reason (like the case you mention). Most of the time, people just want to do some DSP like a low pass filter or resizing.


That’s impressive thanks.

I’m we did Fourier transforms as lab exercise using lasers for my Masters, that was one of quite a few amazing and memorable experiments.

Another was using impedance to pulse lasers, set me to thinking about impedance as a metaphor in software development.

Chemistry A-Level had rates of reaction in Physical chemistry, a few certain key parts being critical in a reaction.


What was the Masters in? Maths or physics or something else? Sounds fun, I wish there'd been more lasers in mine :)


Laser Engineering


To explain further, HT is unable to compose complex valued functions - it is only a complete system for real valued functions, that is: R -> R. So naturally, it's simpler, because it's less expressive.

Let me wow you since HT is not even on the map of exotic or alternative wavelet decompositions.

For example, Hermite transform: https://en.wikipedia.org/wiki/Hermite_transform

or if looking in the discrete space: https://en.wikipedia.org/wiki/Hadamard_transform

See: https://en.wikipedia.org/wiki/Orthogonal_functions

I'm amazed that folks think Fourier series is so special when it's just one of many in a deluge of complete function systems.


the fourier transform is very special indeed because of what its basis functions are the eigenvectors of


So you use one special thing, the Laplacian, to connote the specialty of another? I can go ahead and do the same duality for those other transforms and their respective basis functions. Specialness in mathematics is a fruitless and fealty affair - G.H. Hardy


that's not it


Genuinely curious - How does this save space? If the input was real, the first operation (inside the fft call) makes it complex, which doubles the storage space. Sure the output is real again, but the intermediate storage was 2x the input and output.


Since the input is copied to the real and imaginary components, you don’t actually need to copy it at all.

Think of it like having a pointer (in the C sense) to the real and imaginary image buffers. Set both pointers to the input signal, then run your FFT algo on that.

The key is that it minimizes memory bandwidth, not computation. It’s highly cache coherent.


Thanks, so the real and imag components always map to the same value during the processing, i.e. After every butterfly the real and iamg will be the same value? Cause again if the real and imag components ever differ it means more storage needed.


Correct, they are identical by definition.


What is `area(x)` here? Some kind of numerical integration like the trapezoid rule? Can someone please provide example code which can actually be tested?




Guidelines | FAQ | Lists | API | Security | Legal | Apply to YC | Contact

Search: