DSPRelated.com
Blogs

Add the Hilbert Transformer to Your DSP Toolkit, Part 1

Neil RobertsonNovember 22, 20224 comments

In some previous articles, I made use of the Hilbert transformer, but did not explain its theory in any detail.  In this article, I’ll dig a little deeper into how the Hilbert Transformer works.  Understanding the Hilbert Transformer involves a modest amount of mathematics, but the payoff in useful applications is worth it.

As we’ll learn, a Hilbert Transformer is just a particular type of Finite Impulse Response (FIR) filter.  In Part 1 of this article, I’ll define the Hilbert Transformer, examine the spectrum it creates, and provide some examples of its use.  In part 2, I’ll present a method for designing Hilbert Transformers.

Bandpass Signals

A bandpass signal is one whose spectrum is zero at ω = 0 (Figure 1).  Bandpass signals are typically created by modulating a sinusoidal carrier.  Quadrature down-conversion of such a signal creates a quadrature (I/Q) signal at baseband.  This is a well-known method that is used in demodulation of diverse types of signals.  On the other hand, it is possible to avoid down-conversion and create a bandpass quadrature signal centered at ω0.  It turns out that creating a bandpass quadrature signal is useful for many applications.  The device used to create such a signal is the Hilbert Transformer.

This article is available in PDF format for easy printing

Figure 1.  Spectrum of a Bandpass Signal


Since Hilbert Transformers are frequently used with modulated signals, we start with a definition of bandpass modulation from Sklar [1]:

Bandpass modulation (either analog or digital) is the process by which an information signal is converted to a sinusoidal waveform …  The sinusoid has just three features that can be used to distinguish it from other sinusoids:  amplitude, frequency, and phase.  Thus bandpass modulation can be defined as the process whereby the amplitude, frequency, or phase of an RF carrier, or a combination of them, is varied in accordance with the information to be transmitted.

A modulated signal x(t) can be written as:

$$x(t)=A(t)cos(\omega_0t + \theta(t)) \qquad(1)$$

where A(t) is amplitude (also referred to as the signal envelope), ω0 is carrier frequency, and θ is carrier phase.  An equivalent expression is:

$$x(t)=A(t)cos(\phi(t)) \qquad(2)$$

where $\phi(t)=\omega_0t+\theta(t)$.  We can also represent the signal using a complex exponential:

$$x(t)= A(t)Re\{cos\phi+jsin\phi\} \qquad$$

$$\qquad=Re\{x_r+jx_i\}\qquad\qquad(4)$$

Where $x_r=A(t)cos\phi$  and  $x_i=A(t)sin\phi$.  Figure 2 shows the modulated signal as a vector in the complex plane.

Figure 2.  Modulated Signal in the complex plane.


The Hilbert Transformer

When demodulating x(t), our go-to first step is normally quadrature down-conversion.  This approach allows us to apply lowpass filters to the I and Q channels, which is an important advantage.  Here we’ll take a different approach that does not replace quadrature conversion, but has its own particular uses.  What if we had a device that could generate xr and xi of Equation 4 directly from x(t), without down-conversion?  This would be an alternative first step in demodulation (or other processing) of x(t), and it is what the Hilbert Transformer provides.

We’ll start by considering continuous-time signals.  Figure 3 shows what we are trying to do.  We are looking for a filter h(t) with:

$$input= x_r(t)=A(t)cos(\omega_0t+\theta(t))\qquad\qquad$$

$$output= x_i(t)=A(t)sin(\omega_0t+\theta(t))\qquad(5)$$

We then form the complex signal:

$$x_c(t)= x_r(t)+jx_i(t)\qquad(6)$$



Figure 3.  Hilbert Transformer


Suppose xr is an unmodulated cosine with A(t)= 1.  The Fourier Transform of a cosine is shown in the top of Figure 4.  This spectrum is given by [2]:

$$X_r(\omega)=\pi[\delta(\omega+\omega_0)+\delta(\omega-\omega_0)]\qquad(7)$$

Our goal is to obtain the sine spectrum shown in the bottom of Figure 4 for any value of ω0:

$$X_i(\omega)=H(\omega)X_r(\omega)=j\pi[\delta(\omega+\omega_0)-\delta(\omega-\omega_0)]\qquad(8)$$

where the transfer function H(ω) is the Fourier Transform of h(t).  The H(ω) needed to create Xi(ω) is simply:

$$H(\omega)=\frac{X_i(\omega)}{X_r(\omega)}=\quad j,\;\omega<0$$

$$\qquad\qquad\qquad\qquad\qquad\qquad\qquad -j,\;\omega\ge0\qquad\qquad(9)$$

This rather peculiar-looking transfer function is shown in the middle of Figure 4.  We see that H(ω) provides a +90 degree shift for ω < 0, and a -90 degree shift for ω > 0 (recall that j = ejπ/2).

Once we know the impulse response h(t) that corresponds to H(ω), we can use it to find xi for any modulated bandpass signal.  We can then compute A(t) and φ(t) with reference to Figure 2:

$$ A(t)=\sqrt{x_r^2+x_i^2}\qquad(10)$$

  

$$ \phi(t)=tan^{-1}(x_i/x_r)\qquad(11)$$


Equation 10 represents amplitude demodulation, while Equation 11 represents phase demodulation.  The time derivative of Equation 11 provides dφ/dt, i.e., frequency demodulation.

Figure 4.  Top: Cosine spectrum Xr(ω).  Middle: Hilbert Transformer Transfer Function H(ω).  Bottom: Sine Spectrum Xi(ω) = Xr(ω)* H(ω).


Hilbert Transformer Impulse Response and Structure

The middle plot of Figure 4 shows the frequency response of an ideal continuous-time Hilbert Transformer.  Figure 5 shows the frequency response of an ideal discrete-time Hilbert transformer, where fs is the sample frequency.  We can find the impulse response using the inverse Fourier Transform [3]:

$$h(n)=\frac{2}{n\pi}sin^2(n\pi /2),\quad -\infty < n < \infty \qquad(12)$$

This impulse response is plotted in Figure 6.  We make the following observations:  First, h(n) is odd-symmetric, which follows from the fact that H(f) is pure imaginary.  Second, the sample values are proportional to 1/n, and every other sample of h(n) is zero.  Third, h(n) has infinite length, and it begins before n = 0 (it is non-causal).


Figure 5.  Ideal Discrete-time Hilbert Transformer Frequency Response H(f).

Figure 6.  Ideal Discrete-time Hilbert Transformer Impulse Response.


So how do we design a practical Hilbert Transformer?  There are a few approaches.  One way is to truncate the ideal impulse response and then apply a window function.  Another approach, which we’ll discuss in Part 2, uses a half-band filter as the starting point of the design.  Suppose we use one of these approaches to design a 7-tap Hilbert Transformer.  The coefficients can then be written as:

b = [b-3 0 b-1 0 b1 0 b3]

where, due to the odd symmetry, b1 = -b-1 and b3 = -b-3.  We note that the number of taps is odd, and every other tap value is zero.  The delay of this filter is (ntaps-1)/2 = 3 samples, as opposed to a delay of zero for the non-causal filter.  Because of the non-zero delay, we cannot directly implement the block diagram of Figure 3, in which the xr output is identical to the input.  The Hilbert Transformer can instead be implemented using the FIR structure in Figure 7, where we take xr from the output of the center tap to match the 3-sample delay of xI.  We have also labeled the input x (no subscript), to distinguish it from xr.


Figure 7.  Discrete-time Hilbert Transformer with ntaps = 7.


Hilbert Transformer Output Spectrum

In Equation 6, we defined the Hilbert Transformer’s output xc as complex, which has an effect on its output spectrum.  To examine the spectrum, consider the special case where continuous-time xr(t) is an even-symmetry modulated signal, with carrier frequency ω0.  Then the spectrum Xr(ω) is real and even, as shown in the top of Figure 8.  Applying this signal to the Hilbert Transformer, we get the pure imaginary spectrum Xi(ω), as shown in the middle plot.  We then multiply by j, to get jXi in the bottom plot, which is pure real (recall j2 = -1).  Forming Xc = Xr + jXi, the components centered at - ω0 cancel, and we get the spectrum shown on the right.  We call a signal with a one-sided spectrum an analytic signal.  Note also that the level of the analytic signal’s spectrum is twice that of each component of the real signal’s spectrum.

Now we’ll show that the output spectrum is analytic for a general input signal.  The frequency response of the Hilbert Transformer’s complex output is:

$$G(\omega)=\frac{X_r+jX_i}{X_r}\qquad(13)$$

or,

$$G(\omega)=1+jH(\omega)\qquad(14)$$

Then, from Equation 9,

$$G(\omega)= 0,\quad \omega < 0$$

$$\qquad\qquad\qquad\quad 2,\quad \omega \ge 0 \qquad (15) $$

G(ω), which is pure real, is plotted in Figure 9.  The response is zero for ω < 0, and thus the output is analytic for any input signal.  We again see a factor of two gain for ω > 0.  Note that while the transfer function G(ω) is pure real, the output signal xc(t) is complex.


Figure 8.  Spectrum of Xc(ω) for an even-symmetry bandpass input signal.



Figure 9.  Frequency Response of Hilbert Transformer complex output.


Example 1.  Sanity Check

If the Hilbert Transformer really works, we should be able to use a Matlab model of it to create quadrature sinusoids (cosine/sine) from a sinusoidal input of arbitrary frequency.  But first, we’ll compute the Hilbert Transformer’s frequency response.  We note that H(z) of the Hilbert Transformer is defined with respect to the center tap, so:

$$H(z)=\frac{H_{FIR}(z)}{z^{-K}}\qquad (16)$$

where HFIR(z) is the causal transfer function of the Hilbert Transformer and K samples is the delay of the center tap.

For our model, we’ll use the Hilbert Transformer coefficients listed below (I will discuss Hilbert Transformer design in Part 2).  The length of the filter is 19 taps.  H(f) is calculated from the coefficients b_hilb and b_ct as shown, where b_ct is the impulse response of the center tap and fs = 200 Hz is the sample frequency.

    b_hilb= [-4 0 -21 0 -64 0 -170 0 -634 0 634 0 170 0 64 0 21 0 4]/1024;
    b_ct= [0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0];
    N= 1024;                % number of samples/fft length
    fs= 200;                % Hz sample frequency
    % HT freq response
    [H,f]= freqz(b_hilb,b_ct,N,fs,'whole');

H(f) is plotted over 0 to fs in Figure 10, along with the ideal response.  Note that the frequency range fs/2 to fs (100 – 200 Hz) of the discrete-frequency response corresponds to f < 0 of the continuous response.  We see that the transfer function of a practical Hilbert Transformer departs from the ideal near 0 and fs/2.  There is also some ripple in the passband.  In this case, the response is reasonably accurate over f = 19 to 81 Hz.  Although the transfer function is not ideal, it is nevertheless pure imaginary.  This follows from the odd symmetry of the coefficients.

Figure 10.  19-tap Hilbert Transformer Frequency Response, with ideal response (dashed), fs = 200 Hz.

Now we’ll generate a cosine at 20 Hz, and filter it with the Hilbert Transformer, using the code listed below.  One cycle of the outputs xr and xi are plotted in Figure 11, where you can see the phase difference of 90 degrees.

    fs= 200;                % Hz sample frequency
    Ts= 1/fs;               % s sample time
    f0= 20;                 % Hz frequency of input cosine
    N= 64;                  %    number of samples
    n= 0:N-1;               % time index
    x= cos(2*pi*f0*n*Ts);   % input cosine
    % Filter x using Hilbert Transformer
    x_r= conv(x,b_ct);      % x_r is output of center tap
    x_i= conv(x,b_hilb);
    plot([x_r' x_i'],'.-','markersize',9),grid
    axis([30 40 -1 1])



Figure 11.  One cycle of Hilbert Transformer outputs for input sinusoid at 20 Hz and fs = 200 Hz.

We can also compare the spectra of the input and complex output of the Hilbert Transformer.  We’ll use the same code as above to generate xr and xi, but with changes to the input cosine’s frequency and number of samples:

                    f0= 27;            % Hz frequency of input cosine
                    N= 512;            %    number of samples

Next, we form the complex output xc.  Then we apply a window function [4] to x and xc before performing an FFT:

    x_c= x_r + j*x_i;                % complex output signal
    x_c= x_c(10:end-9);              % limit length to N
    % window x and x_c for fft
    win= kaiser(N,8)'/((N+1)/2);     % Kaiser window function
    x_win= win.*x;
    x_c_win= win.*x_c;
    % spectra of real input and complex output
    X= fft(x_win,N);
    Xc= fft(x_c_win,N);
    X_dB= 20*log10(abs(X));
    Xc_dB= 20*log10(abs(Xc));
    f= n*fs/N;                      % Hz discrete frequency

The resulting spectra are shown in Figure 12.  The output spectrum is approximately analytic.  There is a remaining image component at 173 Hz that is about 60 dB below the desired component.  We also note that, as expected, the magnitude of the analytic signal xc at 27 Hz is 6 dB higher than the input signal components (factor of 2).

Figure 12.  Top: Spectrum of Hilbert Transformer input, f0 = 27 Hz and fs = 200 Hz.

Bottom: Spectrum of Hilbert Transformer complex output xc


Example 2.  Calculating Image Attenuation of a Hilbert Transformer

In figure 12, we see that the Hilbert Transformer image at 173 Hz is attenuated about 60 dB relative to the desired signal at 27 Hz.  Can we find the image attenuation vs. frequency over the entire band from fs/2 to fs?  Indeed we can:  it can be read from a plot of the frequency response of the Hilbert Transformer’s complex output.  From Equation 14, this is:

$$G(f)= 1+jH(f)$$

Where f is discrete frequency.  We can compute H(f) from our Hilbert Transformer coefficients, and then use H(f) to find G(f).

The appendix lists a Matlab function called hilbert_image for calculating and plotting the response of a Hilbert Transformer’s complex output.  For the filter of Example 1, we plot the response using hilbert_image as follows:

    fs= 200;              % Hz sample frequency
    nfft= 512;            %    number of points in freq. response
    b_hilb= [-4 0 -21 0 -64 0 -170 0 -634 0 634 0 170 0 64 0 21 0 4]/1024;
    hilbert_image(b_hilb,nfft,fs);

The Matlab function’s output is shown in Figure 13, where H(f) (top) and G(f) (middle) are plotted using linear amplitude scales.  The response relative to the passband of G is Grel = G/2, which is plotted using a dB amplitude scale (bottom).  Again, the frequency range fs/2 to fs (100 – 200 Hz) of the discrete-frequency response corresponds to f < 0 of the continuous response.  We see from the Grel plot that the response at 173 Hz is about -60 dB, as expected.

Figure 13.  Top: Hilbert Transformer Frequency Response imag[H(f)], fs = 200 Hz.

Middle: Frequency Response of Complex Output, G(f).

Bottom: Frequency Response of Complex Output, Grel = G(f)/2, in dB.


Example 3.  Amplitude Demodulation

An amplitude-modulated signal can be written as:

$$x(t)= M[1+m(t)]cos(\omega_0t)= A(t)cos(\omega_0t)\qquad(17)$$

Where ω0 is the carrier frequency, M is amplitude of the unmodulated carrier, and m(t) is the modulating signal.  If we apply x(t) to a Hilbert Transformer, we can find the modulation envelope A(t) by using Equation 10:

$$ A(t)=\sqrt{x_r^2+x_i^2}$$

The block diagram of a Hilbert AM envelope detector is shown in Figure 14.  The (analog or digital) bandpass filter is needed if noise or adjacent signals are present.


Figure 14.  Envelope detector using the Hilbert Transformer.

Now we’ll create a Matlab model of an AM demodulator.  We’ll use the same Hilbert Transformer coefficients used in the previous examples:

    b_hilb= [-4 0 -21 0 -64 0 -170 0 -634 0 634 0 170 0 64 0 21 0 4]/1024;
    b_ct= [0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0];

Next we’ll create an AM signal with sinusoidal modulation.  Let carrier frequency be 22 Hz, modulation frequency be 3 Hz, and sample frequency be 200 Hz.

    N= 512;                      % number of samples
    fs= 200;                     % Hz sample frequency
    n= 0:N-1;                    % sample index
    Ts= 1/fs;                    % s sample time
    fc= 22;                      % Hz carrier frequency
    fm= 3;                       % Hz modulating frequency
    carr= cos(2*pi*fc*n*Ts);     % carrier
    m= .9*cos(2*pi*fm*n*Ts);     % modulating signal
    A= 1 + m;                    % envelope
    x= A.*carr;                  % modulated signal

To demodulate x, we perform the operations shown in Figure 14:

    % Filter x using Hilbert Transformer
    x_r= conv(x,b_ct);
    x_i= conv(x,b_hilb);
    % demodulation
    y= sqrt(x_r.^2 + x_i.^2);     % y ~= A

The dc restorer is a first-order high-pass IIR filter with coefficients as shown [5]:

    b= [1 -1];
    a= [1 -31/32];
    m_rx= filter(b,a,y);    % m_rx ~= m

The modulated signal x and the demodulated envelope y are plotted in the top of Figure 2, and the recovered signal m_rx is shown in the bottom.  Finally, there are several numerical methods to compute $\sqrt{x_r^2+x_i^2}$ .  An efficient method of moderate accuracy is presented in [6].  Another option is the CoRDiC algorithm [7]. 

Figure 15.  AM Demodulator.  Top: Modulated signal and detected envelope y(t).

Bottom: DC-restored output.


Example 4.  Phase or Frequency Shifter

I described Hilbert-based phase and frequency shifters in an earlier post [8].  The Hilbert Transformer, with its analytic spectrum, is particularly suited here.  Note:  the Hilbert Transformer design in this post uses the window method.  In Part 2 of this article, I’ll describe a different design method that uses half-band filter coefficients as the starting point of the design.


Example 5.  Phase Detector for a Phase-Locked Loop

This is another application described in an earlier post [9].



For the Appendix, see the PDF version of this article.


References

1.  Sklar, Bernard, Digital Communications, Prentice Hall, 1988, p 127.

2.  Hayt, William H., et. al., Engineering Circuit Analysis, 8th Ed., McGraw Hill, 2012, p 765-766.

3.  Lyons, Richard G., Understanding Digital Signal Processing, 3rd Ed., Pearson, 2011, p 487-488.

4.  Robertson, Neil, “The Discrete Fourier Transform and the Need for Window Functions”, DSPRelated.com, Nov. 15, 2021, https://www.dsprelated.com/showarticle/1433.php

5.  Lyons, section 13.23.

6.  Lyons, section 13.2.

7.  Rice, Michael, Digital Communications, a Discrete-Time Approach, Pearson, 2009, section 9.4.

8.  Robertson, Neil, “Phase or Frequency Shifter Using a Hilbert Transformer”, DSPRelated.com, March 25, 2018, https://www.dsprelated.com/showarticle/1147.php

9.  Robertson, Neil, “Digital PLL’s Part 3 – Phase Lock an NCO to an External Clock”, DSPRelated.com, May 27, 2018, https://www.dsprelated.com/showarticle/1177.php



Neil Robertson       November, 2022





[ - ]
Comment by josefseppJanuary 14, 2023

Dear Mr. Niel,

in your article 'Decimator image response' you show how to get this property. Please explain to me what do you mean by 'image response'.

I am following your dsprelated contributions and your publications with interest and I am enthusiastic about your way of explaining the facts.

Wishing you all the best and continued success

Hoffman

[ - ]
Comment by neiroberJanuary 14, 2023

Professor Hoffmann,

I'm glad you like my posts.

As it happens, I deleted that particular post several years ago.  The reason I deleted it was exactly because of the phrase "image response", which was an incorrect usage.  What I was calculating in the post would more correctly be called "alias response".  I would define the alias response as the response to a high-pass filtered impulse, where the HPF is a "brick-wall" filter with cutoff frequency of half the output sample rate.

This response then shows the unwanted aliased spectrum/impulse response after downsampling.

I think the post was mainly correct, except for the usage of the term "image response".

regards,

Neil


[ - ]
Comment by kazJanuary 25, 2023

Thanks Neil for your articles.

decimating filter alias response sounds better though I am used to call it aliasing profile.

I normally check it by creating an imaginary filter same as decimating filter but frequency shifted to Fs/Decimation factor (I believe).

Kaz

%%%%%%%%%%%%%%%%%%%%%%%%%%

%decimate by 2 as example:

h = fir1(30,.6);

%hw shows aliasing profile in slow output domain, after decimation

hw = h.*exp(j*2*pi*(0:length(h)-1)*1/2);

[a1,f]=freqz(h,1,[0:.01:1-.01]*1,1);

[a2,f]=freqz(hw,1,[0:.01:1-.01]*1,1);

plot(f,20*log10(abs(a1)));hold;

plot(f,20*log10(abs(a2)),'r--');

%%%%%%%%%%%%%%%%%%%%%%%%%%%

red line shows aliasing profile

half_band_70698.png

[ - ]
Comment by neiroberJanuary 25, 2023

Hi Kaz,

Very interesting.  Thanks!

Neil


To post reply to a comment, click on the 'reply' button attached to each comment. To post a new comment (not a reply to a comment) check out the 'Write a Comment' tab at the top of the comments.

Please login (on the right) if you already have an account on this platform.

Otherwise, please use this form to register (free) an join one of the largest online community for Electrical/Embedded/DSP/FPGA/ML engineers: