| Hankel (Fourier-Bessel) Transform [message #25846] |
Thu, 19 July 2001 10:53  |
Georg.Pabst
Messages: 2
Registered: July 2001
|
Junior Member |
|
|
Hi,
I'm looking for the Hankel (Fourier-Bessel) Transform, i.e.,
int(0,infinity) f(t)*BesselJ(t*r)t dt being implemented in IDL.
There is a paper "Siegman A. 1980. Quasi fast Hankel transform. Opt.
Lett. 1, 13-15" and one can also find the code in Fortran or C...
Thanks,
Georg
|
|
|
|
|
|
|
|
| Re: Hankel (Fourier-Bessel) Transform [message #25914 is a reply to message #25846] |
Tue, 24 July 2001 12:15  |
thompson
Messages: 584
Registered: August 1991
|
Senior Member |
|
|
Georg.Pabst@nrc.ca (Georg Pabst) writes:
> Hi,
> I'm looking for the Hankel (Fourier-Bessel) Transform, i.e.,
> int(0,infinity) f(t)*BesselJ(t*r)t dt being implemented in IDL.
> There is a paper "Siegman A. 1980. Quasi fast Hankel transform. Opt.
> Lett. 1, 13-15" and one can also find the code in Fortran or C...
> Thanks,
> Georg
Here's an old program that I think might be what you need.
Bill Thompson
FUNCTION HANKEL,F
;
; This function returns the Hankel transform of the argument.
;
S = SIZE(F)
IF S(0) NE 1 THEN BEGIN
PRINT,'*** Variable must be a one-dimensional array, name= F, routine HANKEL.'
RETURN,F
ENDIF
;
X = INDGEN(F)
K = ( 2. * !PI / FLOAT(N_ELEMENTS(X)) ) * X
SC = 0.*X + 1.
IF N_ELEMENTS(SC) GT 3 THEN BEGIN
SC(0) = 3.D0 / 8.D0
SC(1) = 7.D0 / 6.D0
SC(2) = 23.D0 / 24.D0
ENDIF
;
H = BES0( K # X ) # ( K * F * SC )
;
RETURN,H
END
|
|
|
|
| Re: Hankel (Fourier-Bessel) Transform [message #88180 is a reply to message #25914] |
Wed, 26 March 2014 09:25 |
Mats Löfdahl
Messages: 263
Registered: January 2012
|
Senior Member |
|
|
Den tisdagen den 24:e juli 2001 kl. 21:03:48 UTC+2 skrev William Thompson:
> Georg.Pabst@nrc.ca (Georg Pabst) writes:
>
>> Hi,
>
>> I'm looking for the Hankel (Fourier-Bessel) Transform, i.e.,
>> int(0,infinity) f(t)*BesselJ(t*r)t dt being implemented in IDL.
>
>> There is a paper "Siegman A. 1980. Quasi fast Hankel transform. Opt.
>> Lett. 1, 13-15" and one can also find the code in Fortran or C...
>
>> Thanks,
>> Georg
>
> Here's an old program that I think might be what you need.
>
> Bill Thompson
>
>
>
> FUNCTION HANKEL,F
> ;
> ; This function returns the Hankel transform of the argument.
> ;
> S = SIZE(F)
> IF S(0) NE 1 THEN BEGIN
> PRINT,'*** Variable must be a one-dimensional array, name= F, routine HANKEL.'
> RETURN,F
> ENDIF
> ;
> X = INDGEN(F)
> K = ( 2. * !PI / FLOAT(N_ELEMENTS(X)) ) * X
> SC = 0.*X + 1.
> IF N_ELEMENTS(SC) GT 3 THEN BEGIN
> SC(0) = 3.D0 / 8.D0
> SC(1) = 7.D0 / 6.D0
> SC(2) = 23.D0 / 24.D0
> ENDIF
> ;
> H = BES0( K # X ) # ( K * F * SC )
> ;
> RETURN,H
> END
So this thread is from 2001 but it is all I found by googling for "Hankel transform IDL"...
I tried to use the HANKEL function given above but I couldn't make it work. The X=INDGEN(F) makes no sense to me. Should it be X=INDGEN(n_elements(F))? I tried that and as I can't find the BES0 function, I substituted beselj(K # X, 0). This gave me some output, but it completely failed my test of using the same function to compute the inverse Hankel transform and thus getting the original function back. So maybe my changes were all wrong.
So, does anyone know of a useful implementation for IDL? Or C? I just want to convolve some circularly symmetrical functions, but I want to do it as part of a model fitting problem so it would be great if I could save some time by not doing 2D FFTs.
In a response to the post above, Craig commented that it can be done with the discrete Fourier transform. That sounds easy. Is it?
|
|
|
|
| Re: Hankel (Fourier-Bessel) Transform [message #88202 is a reply to message #88180] |
Thu, 27 March 2014 08:12 |
Mats Löfdahl
Messages: 263
Registered: January 2012
|
Senior Member |
|
|
Den onsdagen den 26:e mars 2014 kl. 17:25:06 UTC+1 skrev Mats Löfdahl:
> Den tisdagen den 24:e juli 2001 kl. 21:03:48 UTC+2 skrev William Thompson:
>
>> Georg.Pabst@nrc.ca (Georg Pabst) writes:
>
>>
>
>>> Hi,
>
>>
>
>>> I'm looking for the Hankel (Fourier-Bessel) Transform, i.e.,
>
>>> int(0,infinity) f(t)*BesselJ(t*r)t dt being implemented in IDL.
>
>>
>
>>> There is a paper "Siegman A. 1980. Quasi fast Hankel transform. Opt.
>
>>> Lett. 1, 13-15" and one can also find the code in Fortran or C...
>
>>
>
>>> Thanks,
>
>>> Georg
>
>>
>
>> Here's an old program that I think might be what you need.
>
>>
>
>> Bill Thompson
>
>>
>
>>
>
>>
>
>> FUNCTION HANKEL,F
>
>> ;
>
>> ; This function returns the Hankel transform of the argument.
>
>> ;
>
>> S = SIZE(F)
>
>> IF S(0) NE 1 THEN BEGIN
>
>> PRINT,'*** Variable must be a one-dimensional array, name= F, routine HANKEL.'
>
>> RETURN,F
>
>> ENDIF
>
>> ;
>
>> X = INDGEN(F)
>
>> K = ( 2. * !PI / FLOAT(N_ELEMENTS(X)) ) * X
>
>> SC = 0.*X + 1.
>
>> IF N_ELEMENTS(SC) GT 3 THEN BEGIN
>
>> SC(0) = 3.D0 / 8.D0
>
>> SC(1) = 7.D0 / 6.D0
>
>> SC(2) = 23.D0 / 24.D0
>
>> ENDIF
>
>> ;
>
>> H = BES0( K # X ) # ( K * F * SC )
>
>> ;
>
>> RETURN,H
>
>> END
>
>
>
> So this thread is from 2001 but it is all I found by googling for "Hankel transform IDL"...
I've spent some time trying to understand various matlab codes for calculating Hankel transforms but they seem to require the function to be sampled at some exponential sampling. I need a code that you can call like the fft(), with equidistant sampling.
The OP was happy with porting a matlab code found at www.nmt.edu/~borchers/hankel.html, but that is a code that requires the name of a function that it then evaluates at points of its own choosing. It does not solve my problem.
> I tried to use the HANKEL function given above but I couldn't make it work. The X=INDGEN(F) makes no sense to me. Should it be X=INDGEN(n_elements(F))?
To be specific, X=INDGEN(F), where F is the input function. How is that supposed to work? IDL will accept a vector of dimension lengths, but that vector must be shorter than 8 elements and should anyway consist of integers > 0. That seems like rather limiting requirements.
> I tried that and as I can't find the BES0 function, I substituted beselj(K # X, 0). This gave me some output, but it completely failed my test of using the same function to compute the inverse Hankel transform and thus getting the original function back. So maybe my changes were all wrong.
Here is one test that demonstrates that the transform of the transform does not return the original function (not even times some constant - that would be OK of course):
r = findgen(101)/100.
f = r^3
h = hankel(f)
hh = hankel(h)
cgplot, r, hh, color = 'red'
cgplot, r, f, /over, color = 'blue'
But that is with my edits. Here is the edited function:
FUNCTION HANKEL,F
;
; This function returns the Hankel transform of the argument.
;
S = SIZE(F)
IF S(0) NE 1 THEN BEGIN
PRINT,'*** Variable must be a one-dimensional array, name= F, routine HANKEL.'
RETURN,F
ENDIF
;
;X = INDGEN(F)
X = INDGEN(n_elements(F))
K = ( 2. * !PI / FLOAT(N_ELEMENTS(X)) ) * X
SC = 0.*X + 1.
IF N_ELEMENTS(SC) GT 3 THEN BEGIN
SC(0) = 3.D0 / 8.D0
SC(1) = 7.D0 / 6.D0
SC(2) = 23.D0 / 24.D0
ENDIF
;
; H = BES0( K # X ) # ( K * F * SC )
H = beselj( K # X, 0) # ( K * F * SC )
;
RETURN,H
END
> So, does anyone know of a useful implementation for IDL? Or C?
Nothing? How about FORTRAN?
> In a response to the post above, Craig commented that it can be done with the discrete Fourier transform. That sounds easy. Is it?
One of the matlab functions I found does this. But it won't let me choose my own sampling points.
|
|
|
|
|
|
|
|
comp.lang.idl-pvwave archive
Members
Search
Help
Login
Home
