CODEWITHSUNDEEP
pythonnumpytransformfft
Anton Savostyanov
Anton Savostyanov

Reputation: 53

Python: Fast Hankel Transform for 1d array

I'm trying to find any existing implementation for Hankel Transform in Python (actually i'm more into symmetric fourier transform of two 2d radially symmetric functions but it can be easily reduced to hankel transform).

I do know about hankel python module, but it requires lambda function for input whereas I have only 1d-array.

Any thoughts?

Upvotes: 2

Views: 4217

Answers (4)

VojtaK
VojtaK

Reputation: 678

As of 2023 there is scipy implementation of discrete Hankel transform.

https://docs.scipy.org/doc/scipy/reference/generated/scipy.fft.fht.html

I am also interested in computing inverse 2D Fourier of radially symmetric function. I am quite not sure if the worse is sampling error when I try to sample it on the 2D grid and compute inverse 2D Fourier directly or if I compute Hankel transform and interpolate on logarithmically spaced grid.

Upvotes: 0

EdR
EdR

Reputation: 533

A friend of mine was looking for a Hankel transform in Python and after deciding there was nothing available - partly on the basis of the answers to this question - approached me to share mine.

I have now tidied, documented and released the code as PyHank (docs).

This package performs Quasi-discrete Hankel transforms, as requested by the original question (if a few years too late!) but hopefully this answer will direct future users who find this question (like my friend) to the right place.

Upvotes: 5

StevenMurray
StevenMurray

Reputation: 752

I'm the author of hankel. While I would not recommend to use my code in this case (since as you mention, it requires a callable input function, and its purpose is to accurately compute integrals, not to do a DHT), I will say that it is possible.

All you need to do is interpolate your input 1D array. How to do this is up to you, but typically something like the following works pretty well:

from scipy.interpolate import InterpolatedUnivariateSpline as Spline
import numpy as np

x, y = # Import/create data vectors

# Do this if y is both negative and positive
fnc = Spline(x,y, k=1) #I usually choose k=1 in case anything gets extrapolated.

# Otherwise do this
spl = Spline(np.log(x), np.log(y), k=1)
fnc = lambda x : np.exp(spl(np.log(x)))

# Continue as normal with hankel.transform(fnc, kvec)

The big issue with doing this is in choosing the parameters N and h such that the transform is well approximated for all values of k in kvec. If kvec spans a wide dynamic range, then hankel is very inefficient as it uses the same underlying arrays (of length N) for each k in a transform, meaning the most difficult k sets the performance level.

So again, in short I wouldn't recommend hankel, but if you can't find anything else, it will still work ;-)

Upvotes: 7

Grr
Grr

Reputation: 16079

If you read the documentation for Hankel all the way through the Transforms section you will see that to perform the transformation you call hankel.transform(function, array, ret_err=bool) so you just need the function for whatever form of the transformation you require. I believe there is a list of transformation functions in the Wikipedia entry for Hankel Transformation.

Upvotes: 1

Related Questions