[Scilab-users] Hilbert transform query

[Serge Steer]

X = hilbert(data)  computes the hilbert transform of data(:) not a N dimensionnal Hilbert transform.
Please find attached a function which is intended to compute the N dimensionnal Hilbert transform obtained applying
the 1D hilbert transform to all columns then to all rows, ....

This function exploits the ability of fft to compute all fft along a given dimension with a single call.
Serge
  


Le 08/02/2016 17:51, Lester Anderson a écrit :
> So what would the syntax be for doing a column run Hilbert, and row run Hilbert?
>
> It does sound like the suggesting would do what I think it is meant to.
>
> X = hilbert(data) // does compute pretty fast!
>
> At the moment I have used the netCDF code to read in GMT (net CDF)
> data, which is a 601 x 601 matrix (x,y). So in order to do the column
> and row method is it necessary to transpose the matrix, e.g. x = x' ?
>
> Thanks for any pointers
>
> On 5 February 2016 at 19:02, Tim Wescott <tim at wescottdesign.com> wrote:
>> Any time you go from 1D to 2D you suddenly end up with more than one way
>> to do things, so I'm pretty sure that "how would one..." should really
>> be worded "how would YOU...", or perhaps "how would someone in this
>> field...".
>>
>> It sounds like you want to keep things rectilinear, so it may be best to
>> just apply the transform column-by-column and row-by-row.  That SHOULD
>> work, and if you do it as matrix operations it should be pretty fast in
>> Scilab.
>>
>> On Fri, 2016-02-05 at 12:52 +0000, Lester Anderson wrote:
>>> Hi Serge,
>>>
>>> I am working with grid data, so looking for the 2D Hilbert, the
>>> results I think should appear similar to doing a directional
>>> derivative, where X highlights features with N-S trends and E-W for Y.
>>>
>>> Lester
>>>
>>>
>>>
>>> On 5 February 2016 at 12:24, Serge Steer <Serge.Steer at inria.fr> wrote:
>>>> Le 05/02/2016 10:56, Lester Anderson a écrit :
>>>>> Hello
>>>>>
>>>>> A quick query. How would one define the Hilbert transform of a grid
>>>>> for X and Y directions; looking for two solutions Hx and Hy (for the
>>>>> real values).
>>>> Can you explain more precisely what you expect?
>>>>   Do you want to apply Hilbert transform to each column and to each rows or
>>>> to perform a 2D Hilbert transform?
>>>> Serge
>>>>>
>>>>> Thanks
>>>>> Lester
>>>>> _______________________________________________
>>>>> users mailing list
>>>>> users at lists.scilab.org
>>>>> http://lists.scilab.org/mailman/listinfo/users
>>>>>
>>>> _______________________________________________
>>>> users mailing list
>>>> users at lists.scilab.org
>>>> http://lists.scilab.org/mailman/listinfo/users
>>> _______________________________________________
>>> users mailing list
>>> users at lists.scilab.org
>>> http://lists.scilab.org/mailman/listinfo/users
>> --
>>
>> Tim Wescott
>> www.wescottdesign.com
>> Control & Communications systems, circuit & software design.
>> Phone: 503.631.7815
>> Cell:  503.349.8432
>>
>>
>> _______________________________________________
>> users mailing list
>> users at lists.scilab.org
>> http://lists.scilab.org/mailman/listinfo/users
> _______________________________________________
> users mailing list
> users at lists.scilab.org
> http://lists.scilab.org/mailman/listinfo/users

-------------- next part --------------
function x = hilbertND(x)
    // Marple, S.L., "Computing the discrete-time analytic signal via FFT,"
    // IEEE Transactions on Signal Processing, Vol. 47, No.9 (September
    // 1999), pp.2600-2603
    // http://ieeexplore.ieee.org/iel5/78/16975/00782222.pdf?arnumber=782222

    if  type(x)<>1 then
        error(msprintf(gettext("%s: Wrong type for input argument #%d: Array of floating point numbers expected.\n"),"hilbert",1))
    end
    if x==[] then return;end
    if ~isreal(x,0) then
        error(msprintf(gettext("%s: Input argument #%d must be real.\n"),"hilbert",1));
    end
    d=size(x);
    for k=1:size(d,'*')
      n=d(k)
      if k<>1 then 
        x=permute(x,[k,1]);
      end
      dp=size(x)
      x=matrix(x,d(k),-1)
      x = fft(real(x),-1,1);
      h=zeros(x)
      no2 = int(n/2);
      if ((2*no2) == n) then  // n is even
        h([1,no2+1],:) = 1;
        h(2:no2,:) = 2;
      else // n is odd
        h(1,:) = 1;
        h(2:(n+1)/2,:) = 2;
      end
      x = fft(x.*h,1,1);
      x=matrix(x,dp)
      if k<>1 then x=permute(x,[k,1]);end
    end
endfunction