Tetsuo Sasao and André B. Fletcher
Introduction to VLBI Systems
Chapter 4
Lecture Notes for KVN Students
Partly based on Ajou University Lecture Notes(to be further edited)
Version 1. (Unfinished.)
Issued on February 19, 2006.

Very Long Baseline Interferometry
Contents
1 Technologies Which Made VLBI Possible
1.1 Basics of Digital Data Processing . . . . . . . . . . . . . . . .
1.1.1 Analog Processing Versus Digital Processing . . . . . .
1.1.2 Sampling and Clipping . . . . . . . . . . . . . . . . . .
1.1.3 Discrete–Time Random Process . . . . . . . . . . . . .
1.1.4 Stationary Discrete–Time Random Process . . . . . . .
1.1.5 Sampling . . . . . . . . . . . . . . . . . . . . . . . . .
1.1.6 Comb Function . . . . . . . . . . . . . . . . . . . . . .
1.1.7 Fourier Transform of a Comb Function Is a Comb Function . . . . . . . . . . . . . . . . . . . . . . . . . . . .
1.1.8 Spectra of Discrete–Time Processes . . . . . . . . . . .
1.1.9 Spectra of Sampled Data . . . . . . . . . . . . . . . . .
1.1.10 Inverse Relations for Spectra of Sampled Data . . . . .
1.1.11 Sampling Theorem . . . . . . . . . . . . . . . . . . . .
1.1.12 Optimum Sampling Interval . . . . . . . . . . . . . . .
1.1.13 Sampling Function . . . . . . . . . . . . . . . . . . . .
1.1.14 Correlations of Nyquist Sampled Data with Rectangular Passband Spectra . . . . . . . . . . . . . . . . . . .
1.1.15 S/N Ratio of Correlator Output of Sampled Data . . .
1.1.16 Nyquist Theorem and Nyquist Interval . . . . . . . . .
1.1.17 Higher–Order Sampling in VLBI Receiving Systems . .
1.1.18 Clipping (or Quantization) of Analog Data . . . . . . .
1.1.19 Probability Distribution of Clipped Data . . . . . . . .
1.1.20 Cross–Correlation of 1–Bit Quantized Data:
van Vleck Relationship . . . . . . . . . . . . . . . . . .
1.1.21 van Vleck Relationship in Autocorrelation . . . . . . .
1

4
5
5
6
6
9
10
12
13
14
17
19
20
22
25
28
33
36
37
39
42
46
51

1.1.22 Spectra of Clipped Data . . . . . . . . . . . . . . . . . 51
1.1.23 Price’s Theorem . . . . . . . . . . . . . . . . . . . . . . 55
1.1.24 Cross–Correlation of the 2-Bit Quantized Data . . . . . 59
1.1.25 Cross–Correlation Coefficient of the 2–Bit Quantized
Data . . . . . . . . . . . . . . . . . . . . . . . . . . . . 60
1.1.26 Power and Cross–Power Spectra of 2-Bit Quantized Data 63
1.1.27 Dispersion of Digital Correlator Output . . . . . . . . . 66
1.1.28 S/N Ratio of Digital Correlator Output . . . . . . . . . 67
1.1.29 Optimum Parameters v0 and n in the 2–Bit Quantization 69
1.1.30 Effect of Oversampling in S/N Ratio of Clipped Data . 70
1.1.31 Coherence Factor and Sensitivity with Given Bit Rate 73
1.2 Frequency Standard . . . . . . . . . . . . . . . . . . . . . . . . 76
1.2.1 VLBI Requires “Absolute” Frequency Stability . . . . . 76
1.2.2 How to Describe Frequency Stability? . . . . . . . . . . 78
1.2.3 Types of Phase and Frequency Noises . . . . . . . . . . 80
1.2.4 Time Domain Measurements . . . . . . . . . . . . . . . 82
1.2.5 “True Variance” and “Allan Variance” of Fractional
Frequency Deviation . . . . . . . . . . . . . . . . . . . 82
1.2.6 True Variance and Allan Variance through Power Spectrum of Fractional Frequency Deviation . . . . . . . . . 85
1.2.7 Time–Interval Dependence of Allan Variance . . . . . . 87
1.2.8 True Variance and Allan Variance through Autocorrelation of Phase Noise . . . . . . . . . . . . . . . . . . . 89
1.2.9 Coherence Function . . . . . . . . . . . . . . . . . . . . 89
1.2.10 Approximate Estimation of Coherence Time . . . . . . 91
1.2.11 Estimation of Time–Averaged Phase Noise . . . . . . . 92
1.3 Time Synchronization . . . . . . . . . . . . . . . . . . . . . . 94
1.4 Recording System . . . . . . . . . . . . . . . . . . . . . . . . . 95
2 Overview of the VLBI System101
2.1 MK–3 As a Prototype of Modern VLBI Systems . . . . . . . . 101
2.1.1 Dual–Frequency Reception . . . . . . . . . . . . . . . . 102
2.1.2 First Frequency Conversion . . . . . . . . . . . . . . . 102
2.1.3 Transmission to Control Building . . . . . . . . . . . . 103
2.1.4 Intermediate Frequency Distributer . . . . . . . . . . . 103
2.1.5 Baseband Conversion . . . . . . . . . . . . . . . . . . . 104
2.1.6 Formatter . . . . . . . . . . . . . . . . . . . . . . . . . 106
2.1.7 Data Recorder . . . . . . . . . . . . . . . . . . . . . . . 106
2.1.8 Phase and Delay Calibration . . . . . . . . . . . . . . . 106
2.1.9 Hydrogen Maser Frequency Standard . . . . . . . . . . 108
2.1.10 Automated Operation by Control Computer . . . . . . 1082

2.2 Modern VLBI Systems . . . . . . . . . . . . . . . . . . . . . . 109
2.2.1 New Recording and Fiber–Link Systems . . . . . . . . 109
2.2.2 Digital Baseband Converters . . . . . . . . . . . . . . . 113
2.2.3 e–VLBI . . . . . . . . . . . . . . . . . . . . . . . . . . 114
2.2.4 VLBI Standard Interface (VSI) . . . . . . . . . . . . . 116
3 Difficulties in Ground–Based VLBI

120

4 Correlation Processing in VLBI

120

5 Observables of VLBI

120

3

1

Technologies Which Made VLBI Possible

Low Noise Amplifier

Low Noise Amplifier

Local Oscillator

IF Frequency
Distributer
Video Converter

Local Oscillator

Phase and
Delay
Calibrator
H-Maser Frequency
Standard

Phase and
Delay
Calibrator

IF Frequency
Distributer

H-Maser Frequency
Standard

Video Converter

Formatter

Formatter

Data Recorder

Data Recorder

Correlation Processor

Digital Data Processing

Figure 1: Schematic view of three major technologies and digital data processing which realized VLBI.
There were three major technologies which enabled realization of radio
interferometers with very long baselines, exceeding thousands of kilometers.
They were:
1. high–stability frequency standard,
2. high–accuracy time synchronization, and3. high–speed high–density recording, or super–wideband data transmission
in nowadays.
In addition, all signal processings essential to VLBI, such as time marking,recording, delay tracking, fringe stopping, and correlation processing are done
digitally in modern VLBI systems. In this sense, rapid progress in digitaltechnology in the last decades has formed a fundament of VLBI, as illustrated
in Figure 1.

4

Therefore, we will first examine theoretical bases of digital signal processing, to an extent which is necessary to understand principles and rolesof digital circuitries used in VLBI. Then, we will see basic elements of the
major technologies mentioned above.

1.1

Basics of Digital Data Processing

1.1.1

Analog Processing Versus Digital Processing

For VLBI, digital processing is much more suited than analog prcessing, as
evident from following comparison.
Analog delay circuit

Digital delay circuit

A coaxial cable can be used as adelay cable (Figure 2), if delay issmaller than ∼ 1 µsec (∼ 300 m).
For larger delay, other means suchas sonic wave must be used.
Almost impossible to use for an intercontinental baseline.
Unstabledelay value against environmental
change. Larger error with larger delay. Change of frequency characteristis with mechanical connection anddisconnection of cables.

A ring buffer composed of a large
RAM and shift registers can be used
for digital delay tracking (Figure 2).
Quite stable and highly repeatablein variable environmental conditions.
Does not need any special tuning.
Accuracy of operation is determinedalmost solely by accuracy of clock
time pulses. Easy to cover large intercontinental delays.
Delay is tracked only discretely withsome loss of signal power.

Analog correlation processing

Digital correlation processing

Multiplication and averaging withanalog devices. Operational range islimited by device characteristics. Affected by the environment.

Multiplication and averaging withlogic devices and counters. Stable operation and high repeatability.
Wide dynamic range.

Analog spectrometer

Digital spectrometer

Narrow–band analog BPF’s withsquare–law detectors (“filterbank
spectrometer”). Difficult to adjust
gains and frequency characteristics
of BPF’s (frequency channels). Affected by the environment. Frequency resolution is fixed by passbands of BPF’s.

Composed of logic devices, shift
registers, and counters. Correlation ⇒ Fourier transformation (XF–
type) and Fourier transformation ⇒correlation (FX–type). Stable operation. Identical characteristics of frequency channels. Frequency resolution could be variable.

5input dataring buffer
readaddresstransmission cabledelay cables

write
addressdelay
delayedoutput data

Figure 2: Examples of analog delay circuit using delay cables (left) and digitaldelay circuit using a ring buffer (right).
1.1.2

Sampling and Clipping

Two important achievements in the theory of digital data processing werevital for VLBI. They are the sampling theorem by Shannon (1949), and theclipping theorem by van Vleck and Middelton (1966, original work was doneby van Vleck during World War II).
VLBI data are sampled, and clipped (or digitized) with 1–bit or 2–bit
quantization (Figure 3).
In the followings, we will see how the information of the analog data isessentially restored from the sampled and clipped data. We will also consider
some loss of information accompanied with the digital data processing.
1.1.3

Discrete–Time Random Process

Discrete sequence of variables x[1], x[2], x[3], · · ·, x[i], · · · is called the “random sequence”, or the “discrete–time random process”, if x[i] at any i is arandom variable, i.e., may vary from trial to trial (Figure 4). This is a “discrete version” of the random process continuously varying in time (henceforth, “continuous–time random process”), which we saw in Chapter 3.
We introduce following statistical concepts for the discrete–time random
process.
• Expectation η[i] of a discrete–time random process x[i] (i = 1, 2, 3, · · ·,n, · · ·) is defined by an equation:
η[i] = hx[i]i,(1)where h i stands for an ensemble average defined by a joint probability
distribution of random variables x[i] (i = 1, 2, 3, · · ·, n, · · ·).
6

Analog signal
0tSampling
0t

Clipping (or quantization)1t-1
Bit stream

0011111110000000000111
To digital circuit
Figure 3: Analog–to–digital (A/D) conversion through sampling, clipping,and bit representation. This figure shows a case of 1–bit quantization.

7

x[i]
trial 1
trial 2
trial 3
trial 4
trial 5
1

2

3

4

5

6

7

8

9

10 11

i

Figure 4: A discrete–time random process is a sequence x[1], x[2], x[3], · · ·,x[i], · · ·, whose value at any i is a random variable.
• Autocorrelation R[m, n] of the discrete–time random process x[i] isdefined by an equation:
R[m, n] = hx[m] x∗ [n]i,(2)where symbol {∗ } stands for complex conjugate.

A discrete–time random process x[i] is called the “white noise” if its
autocorrelation satisfies
R[m, n] = h| x[m] |2 i δmn ,(3)

where δmn is Kronecker’s delta symbol:

 1

δmn =  0

(m = n)(m 6= n).

(4)

• Cross–correlation of two discrete–time random processes:
x[1], x[2], x[3], · · · , x[n], · · ·y[1], y[2], y[3], · · · , y[n], · · ·is defined by
Rxy [m, n] = hx[m] y ∗ [n]i.

8(5)

1.1.4

Stationary Discrete–Time Random Process

The concept of stationary random process, which we introduced in Chapter
3 for continuous–time random process, can be transfered to the discrete–timerandom process in the following way (see, for example, Papoulis, 1984).
• Stationary process.
A discrete–time random process x[i] is called “stationary” if the expectation
η[i] = hx[i]i = η,(6)
is a constant independent of i, and if the autocorrelation
R[n + m, m] = hx[n + m] x∗ [n]i = R[m],(7)depends on difference m of arguments only.
In particular, the stationary random discrete–time process is called the“white noise”, if we haveR[m] = R[0] δm0 .

(8)• Jointly stationary processes.
Two discrete–time random processes x[i] and y[j] are called “jointly
stationary” if they are both stationary, and if their cross–correlation
Rxy [n + m, n] = hx[n + m] y ∗ [n]i = Rxy [m],(9)depends on difference m of arguments only.
Similarly to the continuous–time process case, we introduce• correlation coefficient of a zero–mean stationary discrete–time process
x[i]:
R[m]
,
(10)
r[m] =
R[0]and• cross–correlation coefficient of zero–mean jointly stationary discrete–time processes x[i] and y[j]:rxy [m] = q

Rxy [m]Rxx [0] Ryy [0]

,

(11)

where, autocovariance is just equal to autocorrelation and cross–covarianceis equal to cross–correlation, since we assumed zero–mean processes (i.e.
expectations are equal to zero).
91.1.5

Sampling

Let us call the “sampling” an action which makes a discrete–time processby periodically picking up values of a certain continuous–time process with acertain interval of time (“sampling interval”). The discrete–time process thus
created is called the “time–sample” of the original continuous–time process.
If a discrete–time random process x[n] is a time–sample of a continuous–time random process x(t) with a sampling interval T , i.e. if
x[n] = x(nT ),(12)then statistical properties of x[n] are determined by statistical properties(i.e., by probability distribution) of x(t).
• Expectation and autocorrelation of a random time–sample.
If we denote expectation and autocorrelation of a continuous–time random process x(t) as η(t) and R(t1 , t2 ), respectively, then expectation
and autocorrelation of a random time–sample x[i] = x(iT ) are givenbyη[n] = η(nT ),(13)
andR[m, n] = R(mT, nT ),(14)

respectively.
Proof :
1. η[n] = hx[n]i = hx(nT )i = η(nT ).

2. R[m, n] = hx[m] x∗ [n]i = hx(mT ) x∗ (nT )i = R(mT, nT ).
• Cross–correlation of random time–samples.
If we denote cross–correlation of continuous–time random processes
x(t) and y(t) as Rxy (t1 , t2 ), then cross–correlation of random time–samples x[i] = x(iT ) and y[i] = y(iT ) are given by
Rxy [m, n] = Rxy (mT, nT ).
Proof :
Rxy [m, n] = hx[m] y ∗ [n]i = hx(mT ) y ∗ (nT )i = Rxy (mT, nT ).

10

(15)

• Stationary random time–sample.
If a continuous–time random process x(t) is a stationary random process with constant expectation hx(t)i = η and autocorrelationhx(t + τ ) x∗ (t)i = R(τ ), then a time–sample x[n] = x(nT ) is a stationary discrete–time random process with expectation:
η[n] = η,(16)R[n + m, n] = R[m] = R(mT ).

(17)

and autocorrelation:

Proof :
1. Expectation η[n] of the time–sample x[n]
η[n] = hx[n]i = hx(nT )i = η,is a constant independent of argument n.
2. Autocorrelation R[n + m, n] of the time–sample x[n]R[n + m, n] = hx[n + m] x∗ [n]i = hx(mT + nT ) x∗ (nT )i= R(mT ) = R[m],depends on difference m of arguments only.
• Jointly–stationary random time–samples.
If x(t) and y(t) are jointly stationary continuous–time random processes with cross–correlation hx(t + τ ) y ∗ (t)i = Rxy (τ ), then their time–samples x[n] = x(nT ) and y[n] = y(nT ) are jointly stationary discrete–time random processes with cross–correlation:
Rxy [n + m, n] = Rxy [m] = Rxy (mT ).

(18)

Proof :
1. Time samples x[n] and y[n] are both stationary discrete–time random processes, as we saw above.
2. Their cross–correlation
Rxy [n + m, n] = hx[n + m] y ∗ [n]i = hx(mT + nT ) y ∗ (nT )i= Rxy (mT ) = Rxy [m],depends on difference m of arguments only.
111.1.6

Comb Function

Infinite number of delta functions arranged with equal intervals along a horizontal axis, shown in Figure 5, is called the “comb function”. Thus a comb
Ffunction (t; T ) with period T is given in terms of delta functions δ(t) by
an equation:

where k is an integer.

F(t; T ) =

∞
X

k=−∞

δ(t − kT ),

(19)

...

...
. . . -2T -T

0 T

2T 3T 4T 5T 6T 7T 8T 9T. . t.

Figure 5: A comb function.
An alternative expression of the comb function is known in a Fourier
series form as shown below.
1. Let us expand the comb function to a Fourier series within a rangeTTof interval T :
− <t≤
2
2
∞
X

k=−∞

δ(t − kT ) =

∞
X

2πn

a n ei T t ,

n=−∞

for −

T
T
<t≤ ,
2
2

(20)

where n is an integer.
Folllowing the standard procedure of the Fourier series expansion, we
2πmcalculate n–th Fourier coefficient an by multiplying e−i T t , where mis an arbitrary integer, to both sides of equation (20), and integrating
T
T
<t≤ .
them through22
2. Then, we haveT

• LHS =∞X

Z2

k=−∞ T−2T

δ(t − kT ) e−i 2πm
t
T

dt =

Z2

− T2

2πm

δ(t)e−i T t dt = 1,

since t can be equal to kT (t = kT ) only when k = 0 within the
T
T
<t≤ ,
range
2
2
12

T

• RHS =

∞
X

an

n=−∞

since

Z2

ei

2π(n−m)tT

dt = am T,

− T2
T

Z2

ei

2π(n−m)
t
T



 0

 T

− T2

3. Equating both sides, we have an =
∞
X

k=−∞

δ(t − kT ) =if m 6= n

dt = 

if m = n.

1, for any n, and hence
T

∞
2πn
1 X
ei T t ,T n=−∞

for −

T
T
<t≤ .
2
2

(21)

T
T
4. Although we derived this equality in a limited range − < t ≤ , it22
actually holds for wider range of t. In fact, functions in the both sidesof equation (21) do not change if we substitute t with t + mT with anarbitrary integer m. This means that they are both periodic functions
with period T . Therefore, equation (21) holds for the whole range of t,i.e. −∞ < t ≤ ∞. Thus, we have a general relation
∞
X

k=−∞

δ(t − kT ) =

∞
2πn
1 X
ei T t ,T n=−∞(22)

which holds for any t, and, therefore,F(t; T ) =

∞
2πn
1 X
ei T t ,T n=−∞(23)is the alternative expression of the comb function.
1.1.7

Fourier Transform of a Comb Function Is a Comb Function

Fourier transform(t; T ) of argument t with
2π
(Figure 6).
period T is a comb function of argument ω with period
T
Proof :
According to the general formula of Fourier transformation, we haveF̃(ω; T ) =F̃(ω; T ) of a comb function

Z∞

−∞

F

(t; T )e

−iωt

dt =

Z∞

−∞

13

F

∞
2πn
1 X
ei T tT n=−∞!

e−iωt dt

Comb function
with period T

...

...
0

t

T
Fourier transformation

Comb function
with period ω0

...

...
0
ω0 =

ω

2π
Τ

Figure 6: Fourier transform of a comb function of t with period T is a comb
2π.
function of ω with period ω0 =
T
∞
1 X
=
T n=−∞Z∞= ω0δ(ω − n ω0 ),

∞
X

e

t
−i(ω− 2πn
T )

−∞

n=−∞

∞
2π X
2πn
dt =
δ ω−
T n=−∞
T





(24)

where we introduced a notation:
ω0 =

2π
,
T

and used the general formula of the delta function:Z∞e−iωt dt = 2πδ(ω).

−∞

The RHS of equation (24) is nothing but a comb function of ω with a period
2π
:
ω0 =
T
∞
ω0

X

n=−∞

δ(ω − n ω0 ) = ω0 (ω; ω0 ).
F

Thus,
F

(t; T ) ⇔ ω0 (ω; ω0 ),F

where a symbol ⇔ implies a Fourier transform pair.
1.1.8

Spectra of Discrete–Time Processes

We introduce following definitions.
14

(25)

• Power spectrum.

A power spectrum SD (ω) of a stationary random discrete–time processx[n] with autocorrelation R[m] is given by a discrete Fourier transformwith an arbitrary parameter T (Papoulis, 1984):
SD (ω) =∞X

R[m] e−imωT .

(26)

m=−∞

• Cross–power spectrum.

A cross–power spectrum SDxy (ω) of jointly stationary discrete–timeprocesses x[n] and y[n] with cross–correlation Rxy [m] is given by adiscrete Fourier transform with an arbitrary parameter T :SDxy (ω) =∞X

Rxy [m] e−imωT .

(27)

m=−∞

These spectra, as defined by discrete Fourier transforms with an arbitrary
parameter T in equations (26) and (27), are periodic functions of ω with
2π
a period. They have the same forms as Fourier series, with FourierTcoefficients R[m] and Rxy [m], respectively. The spectra SD (ω) and SDxy (ω)
are, in general, dependent on the arbitrary parameter T . Later, for particularcases of sampled discrete–time processes (time–samples), we will choose T
to be equal to their sampling intervals. Then, we will be able to establish
a relationship between a spectrum of a time–sample and a spectrum of its
original continuous–time process.
• Inverse relations.

Autocorrelation R[m] and cross–correlation RDxy [m] of jointly stationary random discrete–time processes in equations (26) and (27) are given
through the power spectrum SD (ω) and cross–power spectrum SDxy (ω)
by inverse relations:πTR[m] =
2π
T
Rxy [m] =2π

ZTSD (ω) eimωT dω,

(28)

Z

SDxy (ω) eimωT dω,

(29)

− Tπ
π
T

− Tπ

which are nothing but the formulae for Fourier coefficients in the series
expansion.
15

Proof :We prove the inverse relation for the power spectrum SD (ω) given inequation (28) only, since a proof for the cross–power spectrum SDxy (ω)(equation (29)) is given just in a similar way.
π

∞
X

T
1. If SD (ω) =R[n] e−inωT , then R[m] =2πn=−∞
In fact,
π

T
2π

ZT

ZT

SD (ω) eimωT dω.

− Tπ

π

SD (ω) eimωT dω =

− Tπ

since

T
ei(m−n)ωT dω = R[m],R[n]2π n=−∞
π

π

ZT

T


 0

− Tπ

π

T
2. If R[m] =2π

SD (ω) e

imωT

−T



 2π

ei(m−n)ωT dω =

ZT

ZT

∞
X

if n = m,otherwise.
∞Xdω, then SD (ω) =R[n] e−inωT .

n=−∞

− Tπ

We first prove this statement for a limited range of ω confinedππwithin an interval − < ω ≤ . Inserting first equation to theT
T
RHS of second equation, we have
π
T

∞
X

∞
X
T Z
0
−inωT
R[n] e=SD (ω 0 )ein(ω −ω)T dω 0 .
2π π
n=−∞
n=−∞
−T

Note here that∞X

0ein(ω −ω)T is a comb function given in equationn=−∞(23), since, introducing a notation ω0 =
∞
X

n=−∞

0

ein(ω −ω)T=∞X

i 2πn (ω 0 −ω)

e ω0

n=−∞
∞
X

= ω0

k=−∞
∞
X

2π
, we haveT= ω0 (ω 0 − ω ; ω0 )F

δ(ω 0 − ω − k ω0 )2π
2π=
δ ω0 − ω − k.
T k=−∞
T
16





(30)

Therefore, we obtainπ∞X

R[n] e−inωT

∞
X

=

n=−∞

ZT

k=−∞− πTSD (ω 0 )δ ω 0 − ω − k

2π
T



dω 0

= SD (ω),
since the delta function in the integrand takes non–zero value whenω0 = ω + k

2π,T

and this condition holds only when k = 0, provided that ω isππ
confined within the interval − < ω ≤ .
T
TNow, if we extend the function SD (ω) to a periodic function with aππ
2π, beyond the initially imposed interval − < ω ≤ ,period ofTTT
we have∞X

R[n] e−inωT = SD (ω),

n=−∞

for any range of ω.
1.1.9

Spectra of Sampled Data

Let us consider discrete–time processes x[n] and y[n], which are time–samples
obtained by sampling jointly stationary continuous–time random processes
x(t) and y(t) with a sampling interval T :x[n] = x(nT ), and y[n] = y(nT ).
Let autocorrelation of x[n], and cross–correlation of x[n] and y[n], be R[m],and Rxy [m], respectively. They satisfy
R[m] = R(mT ), and Rxy [m] = Rxy (mT ),in view of equations (17) and (18). If we choose the arbitrary parameter T
in the power spectrum SD (ω) and the cross–power spectrum SDxy (ω) of thediscrete–time processes x[n] and y[n], as defined in equations (26) and (27),to be equal to the sampling interval T , i.e.,
T = T,

17(31)

then SD (ω) and SDxy (ω) are related to power spectrum S(ω) and cross–power
spectrum Sxy (ω) of the original continuous–time processes x(t) and y(t) by
equations:
SD (ω) =

(32)

∞
1 X
Sxy (ω + k ω0 ),T k=−∞

SDxy (ω) =where ω0 =∞1 X
S(ω + k ω0 ),T k=−∞

(33)

2π
.
T

Proof :
We prove equation (32) for the power spectrum SD (ω) only, sincea proof of equation (33) for the cross–power spectrum SDxy (ω) is
given just in a similar way.
According to equations (17), (26), and (31), the power spectrumof the time–sample x[n] = x(nT ) is given by
SD (ω) =∞X

R[n] e

−inωT

=

n=−∞

∞
X

R(nT ) e−inωT ,n=−∞where T is the sampling interval. Describing the autocorrelation
R(τ ) of the continuous–time stationary random process x(t) interms of the power spectrum S(ω) through inverse Fourier transformation:Z∞
1
0
R(τ ) =
S(ω 0 ) eiω τ dω 0 ,2π−∞

we have
∞

Z
∞
1 X
0
SD (ω) =S(ω 0 ) ein(ω −ω)T dω 0
2π n=−∞

=

1
2π

Z∞

−∞

S(ω 0 )0ein(ω −ω)T dω 0 .

n=−∞

−∞

According to equation (30),

∞
X

∞
X

0

ein(ω −ω)T is a comb function:

n=−∞
∞
X

n=−∞

e

in(ω 0 −ω)T

∞
2π X
2π
=
δ ω0 − ω − k.
T k=−∞
T

18



Therefore,∞
1 X
SD (ω) =T k=−∞Z∞−∞

S(ω 0 ) δ(ω 0 − ω − k

∞
2π
1 X
=
S ω+k
T k=−∞T

where ω0 =1.1.10

2π) dω 0T∞
1 X=
S(ω + k ω0 ),T k=−∞

2π.
T

Inverse Relations for Spectra of Sampled Data

The inverse relation for the power spectrum SD (ω) of a discrete–time stationary random process x[n]:πT
R[m] =
2π

ZT

SD (ω) eimωT dω,− Tπ

as given in equation (28), must yield an autocorrelation which satisfies R[m] =R(mT ), if the process x[n] is a time–sample x[n] = x(nT ) of a continuous–time stationary random process x(t).
Proof :
Substituting equation (32) to the inverse relation, we obtain
π
T

Z
∞
2π
1 X
R[m] =
S(ω + k ) eimωT dω
2π k=−∞ π
T
−T

∞
X

π
+k 2π
T
T

1
=
2π k=−∞ π

Z

0

2π

S(ω 0 ) eim(ω −k T )T dω 0

− T +k 2π
T
π
+k 2π
T
T

∞
1 X
=
2π k=−∞ πZ0

S(ω 0 ) ei(mω T −2πkm) dω 0

− T +k 2π
T
π
+k 2π
T
T

∞
1 X
=
2π k=−∞ π

1
=
2π

Z∞

Z

0

S(ω 0 ) eimω T dω 0

− T +k 2π
T0

S(ω 0 ) eimω T dω 0 = R(mT ).

−∞19

Similarly, we confirm that the inverse relation:πT
Rxy [m] =2π

ZTSDxy (ω) eimωT dω,− Tπ

in equation (29), gives a cross–correlation of time–samples x[n] and y[n]:Rxy [m] = Rxy (mT ).
1.1.11

Sampling Theorem

Shannon (1949) gave a beautiful proof of the sampling theorem, which heformulated as follows:“If a function f (t) contains no frequencies higher than B cps, it is completely
determined by giving its ordinates at a series of points spaced 1/2B seconds
apart.”A mathematical proof showing that “this is not only approximately, but
exactly, true” was given as follows.
“Let F (ω) be the spectrum of f (t). Then
1f (t) =
2πZ∞

F (ω) eiωt dω

−∞

2πB
1 Z
F (ω) eiωt dω,
=
2π
−2πB

since F (ω) is assumed zero outside the band B. If we let
t=

n
,
2B

where n is any positive or negative integer, we obtain
n
f
2B




1
=
2π

2πB
Z

n

F (ω) eiω 2B dω.

−2πB

On the left are the values of f (t) at the sampling points. The
integral on the right will be recognized as essentially the n–th coefficient in a Fourier–series expansion of the function F (ω), taking the interval −B to +B as a fundamental period. This meansthat the values of the samples f (n/2B) determine the Fourier
coefficients in the series expansion of F (ω). Thus they determine
20F (ω), since F (ω) is zero for frequencies greater than B, and forlower frequencies F (ω) is determined if its Fourier coefficients are
determined. But F (ω) determines the original function f (t) completely, since a function is determined if its spectrum is known.
Therefore the original samples determine the function f (t) completely.”
Shannon (1949) mentioned that Nyquist had pointed out the fundamental
importance of the time interval 1/2B seconds in connection with telegraphy,
and proposed to call this the “Nyquist interval” corresponding to the band
B.
Nowadays, we formulate the sampling theorem in a slightly wider form
(Figure 7).
Analog continuous-time signal

01/(2B)

Spectrum of analog continuous-time data-(n+1)B -nB

0t

Sampled signal

νnB (n+1)B0
1/(2B)

Figure 7: Sampling theorem.
Sampling Theorem:
All the information in an analog continuous–time signal with a passband
spectrum limited within a frequency range nB ≤ ν < (n + 1)B, where B isa bandwidth and n ≥ 0 is an integer, can be preserved, provided the signal issampled with the Nyquist interval 1/(2B).
Here we assume a real process with even or Hermitian symmetric spectrum with respect to frequency. Thus, “spectrum” here implies positive frequency part of the spectrum. Sampling frequecncy 2B, with Nyquist interval
1/(2B), is called the “Nyquist rate”.

21tThe proof of the above theorem is given by equations (32) and (33):∞1 X
2π
SD (ω) =,S ω+k
T k=−∞T

SDxy (ω) =



∞
1 X
2π
Sxy ω + k,
T k=−∞T

as illustrated in Figure 8.
In fact, equations (32) and (33), and Figure 8, show
• if positive frequency part of the analog continuous–time spectrum S(ω)
is confined within a passband nB ≤ ν < (n + 1)B, where n ≥ 0 is an
integer, B is bandwidth, and ν is frequency, and• if sampling interval T is equal to the Nyquist interval: T = 1/(2B),then the analog continuous–time spectrum S(ω) is completely preserved in
the spectrum SD (ω) of the sampled data (First and second panels fromthe top of Figure 8). Therefore, all the information of the original analog continuous–time signal is preserved in the sampled data. This proves the
sampling theorem.
Note, however, that the spectrum SD (ω) of the sampled data in a range
0 ≤ ν < B is inverted in frequency compared with the original analog
continuous–time spectrum S(ω), if the integer n is odd (second panel fromthe top of Figure 8).
On the other hand,• if the sampling interval T is larger than the Nyquist interval, i.e., T >1/(2B), as shown in the third panel from the top of Figure 8, or, eventhough T = 1/(2B), if the analog continuous–time spectrum is confinedwithin aB ≤ ν < (a + 1)B, where a is not an integer, as shown in thebottom panel of Figure 8,
then foots of spectral components with different n’s in equations (32) and(33) are overlapped with each other (this is called the “aliasing”). Therefore, information of the analog continuous–time spectrum S(ω) is no longer
preserved in the spectrum SD (ω) of the sampled data.
1.1.12

Optimum Sampling Interval

In order to see that the Nyquist interval is the optimum interval for sampling,
let us consider an analog continuous–time spectrum which is confined within
22

Spectrum of continuous-time data
S (ω)

-3/2T
-3B

-3/2T
-3B

-1/T -1/2T-2B-B-1/T -1/2T-2B-B-1/T -1/2T
-B

-3/2T
-3B

-1/T -1/2T
-2B
-B

0
S (ω)

0
S (ω)

0

1/2T
B

1/2T
B

1/2T

S (ω)

0

1/T 3/2T
2B 3B

1/T 3/2T
2B 3B

1/T
B

Spectrum of sampled data
SD (ω)

ν-3/2T -1/T -1/2T-3B -2B-B

SD (ω)

1/2TB

1/T 3/2T
2B 3B. .ν .

-3/2T -1/T -1/2T-3B -2B-B

0
SD (ω)

1/2TB

1/T 3/2T
2B 3B

...

...

ν

1/T 3/2T
2B 3B

0

...

-3/2T -1/T -1/2T
-B

B

1/2T
B

. . ν.

...

ν

0

1/2T

SD (ω)...

ν-3/2T -1/T -1/2T-3B -2B
-B

0

1/2T
B

1/T
B

3/2T

ν

. . ν.

1/T 3/2T
2B 3B

Figure 8: Four cases of relation between spectrum of analog continuous–timedata and spectrum of sampled data given by equations (32) and (33). Top:
analog continuous–time spectrum is confined within a passband 2mB ≤| ν |<(2m + 1)B and sampling interval T is equal to Nyquist interval T = 1/(2B).
Second from the top: analog continuous–time spectrum is confined withina passband (2m + 1)B ≤| ν |< 2(m + 1)B and T = 1/(2B). Third fromthe top: T > 1/(2B). Bottom: analog continuous–time spectrun is confinedwithin a passband with boundaries of non–integer multiples of B, i.e., aB ≤|ν |< (a + 1)B, and T = 1/(2B). Here we adopted notations, ν: frequency, B:
bandwidth of the analog continuous–time spectrum, m ≥ 0: an integer, anda: a non–integer number. We assume a real process, and, therefore, an evenor Hermitian symmetric spectrum with respect to frequency. All informationof the original analog continuous–time data is completely preserved in thesampled data in the first two cases, but a part of the information is lost after
sampling in the last two cases.

23

0

ν

B

Figure 9: A band–limited baseband spectrum confined within a frequency
range 0 ≤ ν < B. Here, positive frequency part only is shown.
a baseband 0 ≤ ν < B. Here, the “baseband”, or otherwise called the“video–band”, implies a frequency band containing 0 Hz (or “DC”, whichmeans “direct current”) as the lowest frequency, such as shown in Figure 9.
This is a particular case of the passband spectrum within nB ≤ ν < (n+1)B
when n = 0. Also, we assume that the bandwidth B here corresponds to an
actual extent of the spectrum, that means the spectrum is non–zero in theinside of the interval 0 ≤ ν < B, but zero in the outside.
Then, we can conceive three cases which are shown in Figure 10.
Spectrum of continuous-time data
Spectrum of sampled data
1 OversamplingT<
2B-1/2T -B0

B 1/2Tν
T=-B=-1/2T

0 B=1/2T

ν
T>

-B
-1/2T

0

B
1/2T

ν

-1/T

-1/2T

1
Nyquist sampling
2B-2/T -3/2T -1/T -1/2T

1
Undersampling
2B

-3/T

-2/T

-1/T
-B

0

1/2T

1/T

0

1/2T 1/T

3/2T

0

1/T

2/T

B

ν

2/T

3/T

Figure 10: Spectra of sampled data in oversampling T < 1/(2B) (top),
Nyquist sampling T = 1/(2B) (middle), and undersampling T > 1/(2B)(bottom).
1. Oversampling
If we sample an analog continuous–time signal with sampling interval
24

ν

ν

smaller than the Nyquist interval T < 1/(2B), as shown in top panel of
Figure 10, we will have larger number of data points per unit durationof time, but information contained is not improved at all, comparedwith the Nyquist sampling (sampling with Nyquist interval) case shownin the middle panel of Figure 10. Such a samping with an interval
T < 1/(2B) is called the “oversampling”.
2. Undersampling
On the contrary, if the sampling interval is larger than the Nyquist
interval T > 1/(2B), a part of information in the original analog
continuous–time signal is lost in the sampled data due to the aliasing, as shown in bottom panel of Figure 10. This case is called the“undersampling”.
3. Nyquist sampling
Therefore, the Nyquist sampling, i.e. sampling with the Nyquist interval T = 1/(2B), is the optimum sampling for an analog continuous–time signal with a band–limited baseband spectrum (middle panel of
Figure 10).
1.1.13Sampling Function

We saw above that autocorrelation R(τ ) and cross–correlation Rxy (τ ) ofjointly stationary continuous–time random processes x(t) and y(t) with band–
limited baseband spectra are completely restored from autocorrelation R[n] =R(nT ) and cross–correlation Rxy [n] = Rxy (nT ) of time samples of the processes x[n] = x(nT ) and y[n] = y(nT ), provided that the sampling interval
T is shorter than or equal to the Nyquist interval, i.e. T ≤ 1/(2B), where B
is bandwidth of the spectra. But how can we functionally express R(τ ) andRxy (τ ) through R[n] and Rxy [n]?
The answer is given by the so–called “second part of the sampling theorem”, which states that they satisfy equations:R(τ ) =

∞
X

R[n]

n=−∞

and
Rxy (τ ) =∞X

Rxy [n]n=−∞

The sinc function here:

SAn (τ ) ≡hsin Tπ (τ − nT )hπ(τ − nT )Th,sin Tπ (τ − nT )π(τ − nT )T

sin Tπ (τ − nT )π(τ − nT )
T

25

i

i

,

(34)
i

.

(35)

(36)

is called the “sampling function”.
Proof :
We prove equation (34) for the autocorrelation only, since equation (35) for the cross–correlation can be derived exactly in thesame way.
According to equation (32), power spectrum S(ω) of a stationary
random continuous–time process x(t) and power spectrum SD (ω)
of its time–sample x[n] = x(nT ) with a sampling interval T are
related to each other by∞
1 X
2π
SD (ω) =.
S ω+k
T k=−∞T

Therefore, if the sampling interval T is shorter than or equal
to the Nyquist interval, T ≤ 1/(2B), we obtain the full analog
continuous–time baseband spectrum S(ω) by multiplying to thespectrum of the time sample SD (ω) a rectangular window function P (ω) satisfying
P (ω) =
 1

− Tπ < ω ≤ Tπ ,


 0

(37)

otherwise,

as illustrated in Figure 11. In fact, we haveS(ω) = T P (ω) SD (ω),(38)

from equation (32).
Continuous-time spectrum S(w)Spectrum of sampled data SD (ω)
P(ω)

-π/T

0

π/T

ω

-2π/T

-π/T

0

π/T

2π/T

ω

Figure 11: Analog Continuous–Time spectrum S(ω) (left panel) is fully restored from spectrum of time samples SD (ω) (right panel) by multiplying a
rectangular window function P (ω) given in equation (37).

26

If we introduce Fourier transform pairs S(ω) ⇔ R(τ ), P (ω) ⇔p(τ ) and SD (ω) ⇔ RD (τ ), equation (38) implies
R(τ ) = T p(τ ) ∗ RD (τ ) = TZ∞−∞p(τ − α) RD (α) dα,(39)

in view of the convolution theorem which we saw in Chapter 3,
where symbol “∗” stands for the operation of convolution.
We saw in Chapter 3 that inverse Fourier transform of the rectangular window function is a sinc function:
π
T

∞



πτ



sin T
1 Z iωτ
1 Z
P (ω)eiωτ dω =.
e dω =p(τ ) =
2π
2π π
πτ
−∞

(40)

−T

On the other hand, inverse Fourier transform of equation (32)yields1
RD (τ ) =2π=Z∞2π iωτ
)e dω
T

2πk
1
R(τ )
e−i T τ = R(τ )δ(τ − nT ),
T
n=−∞
k=−∞

(41)

−∞

1
dω =
2πT

Z∞ X
∞

S(ω + kSD (ω)eiωτ∞X−∞ k=−∞∞
X

where we used the shift theorem and the property of the comb
function given in equation (22):
∞
∞
X
2πk
1 X
e−i T τ =δ(τ − nT ).
T k=−∞n=−∞

Therefore, equation (39) is reduced toR(τ ) = T∞X

Z∞ sinn=−∞−∞

=

∞
X

n=−∞

R(nT )π(τ − α)Tπ(τ − α)

sin





R(α) δ(α − nT ) dαπ(τ − nT )Tπ(τ − nT )T,which proves equation (34), since R(nT ) = R[n].
Now we are in position to answer to an interesting question: why accuracy of delay determination in VLBI can be much superior (i.e. smaller)27

than a sampling interval of digitized voltage signals, from which the delay isdetermined? For example, typical delay accuracy of Mark III VLBI system
has been 0.1 nanosecond (1 × 10−10 sec), while a typical sampling intervalin the Mark III observation has been 125 nanosecond (= 1/(2 × 2 MHz) ).
Details apart, an essential point of the answer is in the sampling theorem:
Nyquist sampled data are capable of determining the delay as accurately as
continuous–time data, from which the sampled data are formed, since they
are equivalent to each other in view of the sampling theorem.
The 2B optimal rate and the sampling function have been independently
discovered by a number of researchers in different countries, besides Shannon(1949). The history even goes back to the 19th Century. Interested readers
could consult with a review paper by Meijering (2002).
1.1.14

Correlations of Nyquist Sampled Data with Rectangular
Passband Spectra

Let us consider continuous–time stationary random processes x(t) and y(t)with rectangular power spectra S(ω) with a passband of bandwidth B: a

2πnB ≤| ω |< 2π(n + 1)B,

S(ω) =  0otherwise,(42)

where n is an integer, and n = 0 corresponds to the particular case of thebaseband spectrum (Figure 12).
S(ω)
a
...
-(n+1)B -nB

...
-B

0

B

nB (n+1)Bν = ω /2π

S(ω)
a
-B

0

B

ν = ω /2π

Figure 12: Rectangular passband (top) and baseband (bottom) power spectra.
If we sample the data with Nyquist interval T = 1/(2B), then we obtainfollowing properties for correlations of the time samples.
28

Autocorrelation:
Autocorrelation R(τ ) of an original continuous–time process x(t) is obtained by an inverse Fourier transformation of the passband power spectrum
S(ω) (top panel of Figure 12):1
R(τ ) =
2π=

Z∞

S(ω) e

iωτ

−∞

a
<
π

2π(n+1)B
Z

1
dω =
<
π

S(ω) eiωτ dω0

a  i2π(n+ 21 )Bτ< eπeiωτ dω =

2πnB

−iπBτ

iπBτ

"

Z∞

#

1
a
e
−e
< ei2π(n+ 2 )Bτ
π
iτ




sin(πBτ )
1
Bτ ,
= 2aB
cos 2π n +
πBτ
2

=

ZπB

0



eiω τ dω 0 


−πB

(43)

which has the familiar “white–fringe” form with the fringe pattern enclosedby the bandwidth pattern, as we saw in Chapter 3.
In the baseband spectrum (bottom panel of Figure 12), we have n = 0,
and the autocorrelation of the continuous–time process has a sinc functionform:
sin(2πBτ )R(τ ) = 2aB
.
(44)
2πBτ
For the correlation coefficient of the continuous–time process:r(τ ) =

R(τ )
,
R(0)

we have,1sin(πBτ )cos 2π n +πBτ2
in a case of the general passband spectrum, and

r(τ ) =r(τ ) =







Bτ ,

(45)

sin(2πBτ )
,
2πBτ

(46)

in the particular case of the baseband spectrum.
Now, if we sample the continuous–time process x(t) with the Nyquist
interval T = 1/(2B), correlation coefficient of the time sample is given byr[m] = r(mT ) =

sin



mπ
2
mπ
2

29



1
cos mπ n +
2






,

(47)

r (τ)

r (τ)

1

1

0.8
0.6
0.5
0.4
0.2
0

0

-0.2
-0.4
-0.5
-0.6
-0.8
-1

-3

-2

-1

0

1

2Bτ

2

-1
-3

3

-2

-1

0

1

2Bτ

2

3

r [m]

-5

-4

-3

-2

-1

0

1

2

3

4

5

m

Figure 13: Correlation coefficient r(τ ) of a continuous–time process withrectangular passband spectrum when n = 2 (top left) and baseband spectrum(top right), as given in equations (45) and (46). When this process is sampledwith the Nyquist interval T = 1/(2B), correlation coefficient r[m] of the timesample has the “white–noise” form which is equal to 1 if m = 0, and equal
to 0 if m 6= 0 (bottom), since r[m] = r(mT ) = 0 for all m except for m = 0,
as shown in the top two panels.
for the passband spectrum, andr[m] = r(mT ) =sin(mπ),mπ(48)

for the baseband spectrum, in particular. Both equations (47) and (48) show
the “white–noise” form of the time sample:
r[m] = δm0 =



 1

 0

if m = 0,
if m 6= 0,(49)

as given in equation (8), where δij is Kronecker’s delta symbol. This showsthat different sample points are not correlated, and therefore independentof each other, in time samples of Nyquist sampled data with rectangularpassband spectra.
Relationship between the correlation coefficient of the original continuous–time data and that of the sampled data is illustrated in Figure 13.
30

Cross–correlation:
If a cross–power spectrum Sxy (ω) of jointly stationary continuous–timerandom processes x(t) and y(t) is real (i.e., has zero phase), and rectangular with bandwidth B, such as shown in Figure 12, the situation is muchthe same with the autocorrelation case discussed above, and their cross–correlation has the same functional form as equation (43) or (44). Therefore,cross–correlation coefficient of their time samples has the “white–noise” form,
proportional to the one given in equation (49).
Let us now consider a little more general case, when amplitude A(ω) of
the cross–power spectrum Sxy (ω) is rectangular, as given in equation (42),but phase is non–zero due to some delay τd between correlated signals inprocesses x(t) and y(t), which may in general contain both the signals anduncorrelated noises, just like in an interferometer problem. In such a case,the cross–power spectrum, which contains the signal contribution only, hasa form:Sxy (ω) = A(ω) e−iωτd ,(50)
as we saw in Chapter 3.
Strictly speaking, actual passband spectra to be sampled in realistic interferometers are IF spectra after the frequency–conversion, and hence theirphase spectra usually do not cross the origin, i.e., phases are non–zero atω = 0, unlike in equation (50), as we discussed in Chapter 3. Nevertheless,
we adopt equation (50) for simplicity, assuming an idealized case of “RFcorrelation”, or a case when the “fringe stopping” is ideally performed sothat the phase crosses the origin, but the phase slope still remains due to animperfect “delay tracking”.
Then, in view of the shift theorem, the cross–correlation Rxy (τ ):1
Rxy (τ ) =2π

Z∞

Sxy (ω) eiωτ dω,−∞

should have a similar form as given in equation (43) or (44), but argument
τ is replaced by τ − τd . Thus, cross–correlation coefficient is given byrxy (τ ) = q=

Rxy (τ )Rxx (0) Ryy (0)






sin[πB (τ − τd )]
1


ρ
cos
2π
n
+
B
(τ
−
τ
)

d


πB (τ − τd )
2






sin[2πB (τ − τd )] ρ

2πB (τ − τd )31

(passband),(51)(baseband),where ρ is the maximum cross–correlation coefficient:
ρ= q

Rxy (τd )
Rxx (0) Ryy (0)

.

(52)

Note that the cross–correlation coefficient in the case of the passband spectrum with n 6= 0 again shows the “white–fringe” form with the cosine “fringe
pattern” enclosed within the sinc function envelope of “bandwidth pattern”.
rxy(τ) / ρrxy[m] / ρ rxy(τ) / ρ

1

1

0.5

0.5

0

0

-0.5

-0.5

rxy[m] / ρ

τd
-1
-3

-2

-1

0

2Bτ

τd
1

2

3

m

-1

-3

-2

-1

0

2Bτ

1

2

3

m

Figure 14: Cross–correlation coefficient rxy (τ ) of jointly stationary
continuous–time processes x(t) and y(t) with rectangular passband spectrum when n = 2 (left: dotted line) and baseband spectrum with n = 0(right: dotted line), as given by equation (51). τd is a time delay between
correlated processes x(t) and y(t). Also shown by bars is cross–correlation
coefficient rxy [m] = rxy (mT ) of time samples x[i] = x(iT ) and y[j] = y(jT )sampled with the Nyquist interval T = 1/(2B), where B is a bandwidth of
the rectangular spectrum (equation (53)). Horizontal axis shows the timedifference τ normalized by the sampling interval T = 1/(2B). Vertical axis
is the cross–correlation
coefficient rxy (τ ) normalized by its maximum value:qρ = Rxy (τd )/ Rxx (0) Ryy (0), where Rxy (τ ), Rxx (τ ), and Ryy (τ ) are cross–correlation and autocorrelations of x(t) and y(t), correspondingly.
Therefore, if we sample x(t) and y(t) with the Nyquist interval T =1/(2B), cross–correlation coefficient of the time samples is
rxy [m] = rxy (mT )

32

h 
i

τd
π
 



sin
m
−

1
τd

2
T




ρ π
cos π n +
m−


2
T

m − τTd


2


= 


h 
i

τd


sin
π
m
−

T




ρ


τd

π m−

(passband),(53)(baseband).

T

Relationship between the cross–correlation coefficient of the original continuous–time data and that of the sampled data is illustrated in Figure 14.
The cross–correlation coefficient rxy [m] of the time samples given by equation (53) no longer has the symmetric “white–noise” form, as shown in equation (49) and in the bottom panel of Figure 13, due to the parallel shift of thecross–correlation coefficient of the continuous–time data along the horizontal
axis which is caused by the delay τd . Also, it now depends upon n, i.e. uponlocation of the passband spectrum on the frequency axis, since the “fringe
pattern” in the “white–fringe” depends on the location.
Thus, in the cross–correlation coefficient, the simple “white–noise” form
and the independence of sample points is obtained only when the delay τd isreduced to zero (τd = 0) by a suitable compensating operation, such as the“delay tracking” in the interferometry.
1.1.15

S/N Ratio of Correlator Output of Sampled Data

Let us now imagine a “semi–analog” correlator (non–existing in reality),
which would multiply and integrate (i.e. time–average) sampled but not
quantized (not clipped) data streams from two antennas of an interferometer.
We will estimate here a signal–to–noise ratio of such a correlator, before
examining actual digital correlators which deal with sampled and quantizeddata.
Let us assume that the two sampled data streams x[i] and y[i] are timesamples of jointly stationary continuous–time random processes x(t) andy(t), which obey the second–order Gaussian probability distribution, as we
assumed in the signal–to–noise–ratio discussion in Chapter 3. We further
assume that x(t) and y(t) have identical rectangular passband spectra withbandwidth B, as given in equation (42), and they are sampled with the
Nyquist interval T = 1/(2B). Then, we have x[i] = x(i T ) and y[i] = y(i T ).
Also, we assume that the delay tracking and the fringe stopping are perfectely
performed beforehand, so that the two input data of exactly the same wave
front are being correlated.
In this case, “correlator output” Rs of the sampled data streams is an33

average of products of time samples over a certain number N :N
1 X
x[i] y[i].
Rs =
N i=1

(54)

Expectation of this correlator output is nothing but the cross–correlation
Rxy [0] of x[i] and y[j] at zero argument, sincehRs i =

N
N
1 X
1 X
hx[i] y[i]i =hx(i T ) y(i T )i = Rxy (0) = Rxy [0].
N i=1
N i=1

(55)

On the other hand, dispersion of this correlator output σs2 is given by
σs2 = hR2s i − hRs i2 ,(56)

as we saw in Chapter 3. hR2s i is described through a double sum of the
fourth statistical momentum in view of equation (54). The fourth statistical momentum is decomposed into a sum of products of second statistical
momenta (correlations), as we discussed in Chapter 3, since x[i] = x(i T )and y[j] = y(j T ) obey the joint Gaussian probability distribution. Thus, wehavehR2s iN X
N
1 X
hx[i] y[i] x[j] y[j]i=
N 2 i=1 j=1

=

N X
N
1 X
{ hx[i] y[i]i hx[j] y[j]i
N 2 i=1 j=1+ hx[i] x[j]i hy[i] y[j]i + hx[i] y[j]i hy[i] x[j]i }

N X
N n
o
1 X
2
=
R
[0]
+
R
[i
−
j]
R
[i
−
j]
+
R
[i
−
j]
R
[j
−
i]
xx
yy
xy
xy
xy
N 2 i=1 j=1
1
1 2= hRs i2 + Rxx [0] Ryy [0] + Rxy
[0]
N
N
1
(57)
= hRs i2 + Rxx [0] Ryy [0] (1 + ρ2 ),N

where ρ = q

Rxy [0]=q

Rxy (0)is the maximum cross–correlation
Rxx [0] Ryy [0]Rxx (0) Ryy (0)coefficient, given in equation (52), in our assumed case with τd = 0. In deriving last two lines of equation (57), we used the “white–noise” relations for
autocorrelations:
Rxx [i − j] = Rxx [0] δij , and Ryy [i − j] = Ryy [0] δij ,
34

(58)

in view of equation (49), and for cross–correlation:
Rxy [i − j] = Rxy [0] δij ,(59)

which is also satisfied since we assumed τd = 0.
Therefore, the dispersion of the correlator output in equation (56) is now
given by1
Rxx [0] Ryy [0] (1 + ρ2 ),
(60)
σs2 =
N
and we obtain the signal to noise ratio SN R:
SN R =√√hRs i
ρ
Rxy [0]
N = √
N.
= q
2
σs
1+ρ
Rxx [0] Ryy [0] (1 + ρ2 )(61)

In the expression of the maximum cross–correlation coefficient ρ, the autocorrelations Rxx (0) and Ryy (0) are usually dominated by system noise contributions from antenna–receiver systems of a radio interferometer, while
the cross–correlation Rxy (0) contains contribution of the signal from a radio source only, as we discussed in Chapter 3. Thus, when we observe acontinuum spectrum source, ρ is approximately given byρ =s

TA 1 TA 2,
TS 1 TS 2(62)

as we saw in Chapter 3, where TA1 , TA2 are antenna temperatures, which
are assumed constant throughout the frequency band B in the case of thecontinuum spectrum source, and TS1 , TS2 are system noise temperatures, ofantenna 1 and antenna 2.
For most of radio sources, TA  TS , and, therefore, ρ  1. In this case,equation (61) is reduced to√
hRs i
SN R ==ρ N =σss

TA 1 TA 2 √N.
TS 1 TS 2(63)

If we denote an integration time of the correlation processing as τa , thenumber of samples N with Nyquist interval 1 / (2 B) is equal to
N = 2 B τa .

(64)

Therefore, equation (63) for the continuum spectrum source is reduced to
SN R =s

TA 1 TA 2 q
2 B τa ,TS 1 TS 2
35(65)which is just identical with what we derived for correlator output of continuous–time voltages in Chapter 3.
This means that the Nyquist sampling does not cause any loss of signal–
to–noise ratio of the correlator output, compared with the continuous–timecase, as expected from the sampling theorem. This also means that there is
no room for the oversampling, with a sampling rate higher than the Nyquist
rate (T < 1 / 2B), in improving the signal–to–noise ratio, despite increased
number of data points. Thus, the Nyquist sampling is really optimum forthe radio interferometry.
√
Note that equation (63) can be interpreted as showing N –fold improvement of the signal–to–noise ratio after repeating and averaging N “measurements” of a power (product of two data streams, in our case). This meansthat measurements of a power made at the Nyquist interval are independent of each other, in the case of the rectangular passband spectra. Thisis a consequence of the independence of time samples themselves discussed
earlier.
1.1.16

Nyquist Theorem and Nyquist Interval

The Nyquist theorem (Nyquist, 1928), which we saw in Chapter 2, says that
thermal noise power Wν per unit bandwidth emitted by a resistor in a thermal
equilibrium with a temperature T is equal toWν = k T,(66)

in the classical limit hν  kT , where k and h are the Boltzmann and rhe
Planck constants, respectively. Therefore, energy E emitted within a rectangular band with a bandwidth B during a time interval t is
E = B t k T.

(67)

1
Since energy per one degree of freedom is equal to k T under the thermal
2equilibrium, number of degrees of freedom in this energy must be NF = 2 B t.
On the other hand, we have NI = 2 B t Nyquist intervals during the timet, for the bandwidth B. In the case of the rectangular band, one Nyquist
interval contains one independent sample, as we saw earlier. Therefore, wehave NI independent samples in the emitted energy during the time t.
The equality NI = NF = 2 B t means that one independent sample(Nyquist interval) in the information theory corresponds to one degree of
freedom in the physics, in the thermal noise.

36

1.1.17

Higher–Order Sampling in VLBI Receiving Systems

In digital data processings as applied to radio astronomy, sampling of received voltage signals has been traditionally done at the basebands (or thevideo–bands), containing DC (zero frequency) as the lowest frequency, after
frequency conversions. This was the safest way for reliable sampling, whenclock rates of sampler circuits were not high enough, and not very stable.

Figure 15: Diagram of receiving system in KVN (Korean VLBI Network)
adopting the higher–order sampling technique (figure brought from KVN
webpage http://www.trao.re.kr/ kvn/).
However, it is not easy, in existing analog filtering technology, to implement a good enough lowpass filter with sharp rectangular edges. Thissituation has often resulted in rather poor frequency characteristics of thebaseband spectra, and made it difficult to achieve high signal–to–noise ratio,
close to the one expected from an ideally rectangular spectrum. Also, because of this difficulty, high quality baseband converters tend to be expensive,
especially when wide frequency bands are required.
37

Recently, Noriyuki Kawaguchi and his colleagues successfully applied so–called “higher–order sampling” technique to a number of VLBI systems,
including Japanese VERA (VLBI Exploration of Radio Astrometry). The
higher–order sampling is the sampling at a passband with n > 0, discussed
earlier. In general, it is easier to design good analog bandpass filters, withnearly rectangular band shapes, when the ratio B/ν0 of bandwidth B to central freuency ν0 is smaller. Therefore, it is easier to make a nearly rectangularwideband filter for a passband, than for a baseband. In fact, the higher–ordersampling technique has been effective in wideband receiving systems withtypical bandwidth of 512 MHz or wider, for realizing better frequency characteristics and higher signal–to–noise–ratio (Iguchi and Kawaguchi, 2002).
Figure 15 shows a diagram of KVN (Korean VLBI Network) receiving
system which adopts the higher–order sampling technique. A “baseband
converter” cuts off a 512 MHz band from a 2 GHz–wide first IF signal with
8.5 GHz center frequency (i.e. 7.5 – 9.5 GHz band), and converts it to a 512
MHz–wide passband signal with 768 MHz central frequency (i.e. 512–1024
MHz band), which is then sampled by a high–speed sampler.
Note here that, for the bandwidth B = 512 MHz, the passband nB ≤ ν <(n + 1)B with nB = 512 MHz means an odd number of n (i.e. n = 1), whereν is the frequency and n is an integer. Therefore, spectrum of sampled signal
is inverted with respect to the spectrum of original continuous–time signal,
as we saw in Figure 8. In order to avoid possible inconveniences with theinverted spectra, LO (local oscillator) frequency of the “baseband converter”
is chosen so that the passband spectrum is obtained in the lower sideband(LSB), i.e. the spectrum is inverted with respect to the first IF spectrum. In
this way, one can obtain a spectrum of the sampled signal which is identicalwith the one contained in the first IF signal.
It seems worthwhile to mention an interesting question here. The sampling theorem says that the optimal sampling rate for a passband nB ≤ ν <(n + 1)B is 2B. Does this mean that we can use an inexpensive low–speed
4 Msps (mega sample per second) sampler for sampling a signal in a high–frequency passband, say, 10.000 – 10.002 GHz? The answer is “NO”, as N.
Kawaguchi clearly explains. Although the required sampling interval is really
1/(2B), sampling timings must be controlled with much greater accuracy,
better than 1/ν0 , where ν0 is the central frequency of the passband (the “carrier” frequency). Otherwise, we will get all chaotically “jittered” data at eachof VLBI stations, from which we will never find any good fringe. Therefore, a
required sampler must be as accurate as, and as stable as, a 20 Gsps sampler,say. Consequently, the high–speed sampling technology is indispensable forsuccessful application of the wideband higher–order sampling technique.
38

1.1.18

Clipping (or Quantization) of Analog Data

The sampling replaces data that are continuous in time with those discretein time. This is undoubtedly a big step towards digitizing analog data.
However, the sampling alone still leaves values of time samples analog, which
may vary arbitrarily from sample to sample. For a complete analog to digital(A/D) conversion, we need to replace each continuously variable value with
an element of a finite set of discrete values expressible by a certain numberof bits. This step is called the “clipping” or the “quantization”.
Number of discrete values expressible by a given number of bits determines “number of levels of quantization”. Therefore, 1–bit, 2–bit, · · ·, n–bit
quantizations usually correspond to 2–level, 4–level, · · ·, 2n –level quantizations, respectively. There have been exceptions of the n–bit, 2n –level law ofquantization, such as 2–bit, 3–level quantization. However, such exceptional
quantization schemes are rarely used in present–day VLBI systems.
We will denote henceforth a discrete–time process with quantized values
as x̂[i]. If number of quantization levels is m, with discrete values x1 , x2 , · · ·,xm , then x̂[i] must take one of these m values as illustrated in Figure 16. A
quantized process x̂[i] is supposed to be related in a prescribed way to an
original discrete–time process x[i] with analog values before clipping.
x[i]
x8
x7
x6
x5
x4
x3
x2
x1

i

Figure 16: An image of a quantized discrete–time process.
The larger the number of quantization levels (i.e. bits), the more information can remain after the clipping. For reducing data size and increasing
processing speed, however, smaller number of bits is preferable. Thus one
has to choose an optimal number of quantization levels for one’s particular purpose. In VLBI, 1–bit (2–level) and 2–bit (4–level) quantizations aremostly used.
Figure 17 and Table 1 show how quantized values are related to original analog values in cases of VLBI 1–bit (left panel)and 2–bit (right panel)39quantizations.
In 1–bit quantization scheme, clipped quantity takes only one of twovalues: +1 and −1, depending upon if original analog value is positive or
negative, respectively. It is generally accepted that bit 0 is assigned to −1state and bit 1 is assigned to +1 state.
x (clipped)

x (clipped)
+n

+1
+1
x

x
-1

-1
-n
-v0

0

bit

0

sign
bit
mag.
bit

1

0

+v0

0 (1)

0 (1)

1 (0)

1 (0)

0 (1)

1 (0)

0 (0)

1 (1)

Figure 17: Relations between analog (unquantized) values x and clipped(quantized) values x̂ for 1–bit (left) and 2–bit (right) quantizations, respectively. Bit assignments are shown in bottom panels. In case of 2–bit quantization (left), a representative bit assingment for data recorder is shown alongwith that for a correlator chip given in parentheses.

Quantization

1–bit (2–level)

2–bit (4–level)

Analog value

x<0

Clipped valuex̂ = −1x̂ = +10≤xx < −v0−v0 ≤ x < 00 ≤ x < v0x̂ = +1x̂ = +n

Recorder bit ass.

01

s 0, m 0

s 0, m 1

s 1, m 0

s 1, m 1

Correlator bit ass.

01

s 1, m 1

s 1, m 0

s 0, m 0

s 0, m 1x̂ = −nx̂ = −1v0 ≤ x

Table 1: Clipping criteria and bit assingments of VLBI 1–bit and 2–bit quantizations. In the bit assignments for the 2–bit quantization scheme, “s” and“m” stand for sign bit and magnitude bit, respectively.
In 2–bit, 4–level quantization scheme, 4 quantization states are separated
by three threshold values: −v0 , 0, and v0 . They correspond to 4 clippedvalues: −n, −1, +1, and +n, as shown in left panel of Figure 17 and Table
1. Values of parameters v0 and n are chosen so that signal–to–noise ratio40

of the correlator output is maximized. Note that n thus determined is notnecessarily an integer.
Bit assignments to the 4 quantization states are not uniquely standardized
yet in the 2–bit quantization scheme. Existing VLBI data recording systems
mostly adopt 00, 01, 10, and 11 assignments for −n, −1, +1, and +n states,
where first and second bits stand for sign and magnitude bits, respectively,
but a widely used 2–bit quantization correlator chip adopts 11, 10, 00, and
01 assignments, as shown in Figure 17 and Table 1. Therefore, recorded bits
should be rearranged before the correlation, when we use a correlator withthe chip.
Sign bits in the 2–bit quantized data are equivalent to 0 and 1 bits in
the 1–bit quantized data. Therefore, it is usually not difficult to cross–correlate 1–bit quantized and 2–bit quantized data, which are obtained indifferent stations in the same VLBI observation with the same sampling
rate, using a 1–bit correlation mode, if we sacrify some of information in the
2–bit quantized data. In this case, the only necessary thing is to pick up sign
bits from the 2–bit qunatized data and cross–correlate them with the 1–bit
quantized data. Of course, direct cross–correlation of data with different
quantization schemes, such as 1–bit and 2–bit, can be performed withoutloosing information (see, for example, Hagen and Farley, 1973), provided
that a special logical circuit for this purpose is built in a correlator.
Analog signal
+v

0

0

t

-v

0

2-bit quantized signal
+n
+1
t

-1
-n

Figure 18: A schematic view of time variation of original analog data (top)and 2–bit quantized data (bottom).
In the 1–bit quantization scheme, only sign information of the original
analog data remains in the clipped data, and no information on magnitude(amplitude) is left at all, as illustrated in Figure 3. In the 2–bit quantized41

data, some information on the magnitude (amplitude) of the analog data isleft, as illustrated in Figure 18, but the information is very vague.
Then, how much can we restore scientific information contained in theoriginal analog data of VLBI after clipping them with 1–bit or 2–bit quantization scheme, which looks quite rough at least at the first glance?
The clipping theorem gives a surprising answer. The theorem was originally developed by J.H. van Vleck in a study of radar–jamming during theWorld War II conducted by the Radio Research Laboratory of Harvard University (Report No.51 on July 21, 1943), and was made public more than 20
years later by van Vleck and Middleton (1966).
1.1.19Probability Distribution of Clipped Data

Before proceeding to the clipping theorem, we will examine how we candescribe probability distribution of clipped data.
As we discussed earlier, signals from astronomical radio sources, as well as
system noises produced in antenna–receiving systems and in environments,
are well approximated by Gaussian random processes. Therefore, let us consider the data as jointly stationary continuous–time random processes x(t)and y(t), which obey the second–order Gaussian probability density, introduced in Chapter 3. Here we assume a zero–mean case (i.e. expectationsof x(t) and y(t) are equal to zero), and use notations suited to the current
discussions, slightly changing those adopted in Chapter 3. Also, we assume
that both x(t) and y(t) are real processes.
Then, we describe the zero–mean second–order Gaussian probability density of jointly stationary continuous–time random processes x(t + τ ) and y(t)as:f (x, y; τ ) ≡ f (x, y; t + τ, t)=12πσx σyq

2 (τ )
1 − rxy

xy
1
x2
y2
−
2
r
(τ
)
−
+
xy
2
2
σx σy σy2 , (68)e 2[1 − rxy (τ )] σx"#

where we introduced notations: σx2 ≡ Rxx (0), σy2 ≡ Ryy (0), for dispersions ofx(t) and y(t), and
Rxy (τ )rxy (τ ) ≡ q,
Rxx (0) Ryy (0)for cross–correlation coefficient of x(t+τ ) and y(t). Here, Rxy (τ ), Rxx (τ ), andRyy (τ ) are cross–correlation and autocorrelations of x(t) and y(t), as before.
We first assume that the cross–correlation coefficient rxy (τ ) is smaller than42

2
1 in absolute value (i.e. rxy(τ ) < 1) in order to avoid possible singularity in2
our calculations which may ocur when rxy(τ ) = 1.

x(t)
x+dx
x

t+τ

y(t)

time

y+dy
y
time

t

Figure 19: In case of continuous–time processes, Joint probability density
describes probability of the processes to take values contained within infinitesimal ranges at certain periods of time.
Of course, this joint probability density satisfies the general definition of
the cross–correlation:
Rxy (τ ) = hx(t + τ ) y(t)i =Z∞ Z∞

x y f (x, y; τ ) dx dy,

(69)

−∞ −∞

as we saw in Chapter 3 (note that in our zero–mean case the cross–covariance
Cxy (τ ) is just equal to the cross–correlation, i.e. Cxy (τ ) = Rxy (τ ) ).
It is worth to recall here that the joint probability means“f (x, y; τ ) dx dy is a probability for x(t + τ ) and y(t) to bewithin ranges:x ≤ x(t + τ ) < x + dx,y≤y(t)< y + dy,
for any t (see Figure 19). The condition “for any t” corresponds
to our case of the stationary random processes.
Now if we consider discrete–time processes x[i] and y[i], which are timesamples of the above jointly stationary real continuous–time random processes x(t) and y(t), i.e. x[i] = x(iT ) and y[i] = y(iT ), where T is a sampling
43

interval, the time samples x[i] and y[i] are also jointly stationary as we saw
earlier, and their cross–correlation:
Rxy [m] = hx[n + m] y[n]i = Rxy (mT ) = hx(nT + mT ) y(nT )iis described by the same joint probability density of the continuous–timeprocesses f (x, y; τ ) as
Rxy [m] = Rxy (mT ) =Z∞ Z∞

x y f (x, y; mT ) dx dy.

(70)

−∞ −∞

Let us then consider that we clip the time samples x[i] and y[i], andobtain clipped discrete–time processes which we denote as x̂[i] and ŷ[i]. They
now take only discrete values of a finite number N (this means N –level
quantization) x1 , x2 , · · ·, xN and y1 , y2 , · · ·, yN .

^
x[i]
x5
x4
x3
x2
x1
n+m

^
y[i]

i

y5
y4
y3
y2
y1
n

i

Figure 20: Joint probability P (xi , yj ; m) of clipped processes x̂[n + m] andŷ[n] describes probability for them to take certain discrete values xi and yj
on quantization levels (here 5–level case is shown).
In this case, their cross–correlation Rx̂ŷ [m] = hx̂[n + m] ŷ[n]i is decribedby an equation:Rx̂ŷ [m] =

N X
N
X

xi yj P (xi , yj ; m),(71)

i=1 j=1

where P (xi , yj ; m) is a joint probability for x̂[n + m] and ŷ[n] to take values:x̂[n + m] = xi ,ŷ[n] = yj ,44

for any n (See Figure 20).
This cross–correlation is what should be yielded by a digital correlator.
In fact, in view of the ergodicity, the cross–correlation is well approximatedby a time–average of products of N quantized data x̂[n + m] and ŷ[n]:Rx̂ŷ [m] ∼=N
1 X(x̂[m + n] ŷ[n]),
N n=1

(72)

if N is sufficiently large; and the digital correlator is nothing but a machine
which time–averages a large number of products of digital (i.e. sampled andclipped) data.
At the first glance, equations (70) and (71) look quite different, and it
appears difficult to relate Rx̂ŷ [m] with Rxy [m]. However, if the quantization
condition is clearly specified, we can calculate P (xi , yj ; m) rather easily fromf (x, y; mT ) (van Vleck and Middleton, 1966).
In case of the 1–bit, 2–level quantization, the condition is


 +1

for x[i] = x(iT ) ≥ 0,x̂[i] = for x[i] = x(iT ) < 0,ŷ[i] = for y[i] = y(iT ) < 0.

 −1


 +1

 −1

for y[i] = y(iT ) ≥ 0,

(73)

Therefore, the probability for x̂[i] to be +1, for example, is nothing but theprobability for x(iT ) to be 0 ≤ x(iT ) < +∞. Thus we can describe joint
probabilities for all combinations of quantization states through equations:
P (+1, +1; m) =

+∞
Z

 +∞
Z
0





Z0  Z00

P (−1, −1; m) =P (+1, −1; m) =P (−1, +1; m) =f (x, y; mT ) dxf (x, y; mT ) dx dy,

−∞ −∞
 +∞
Z
Z0




−∞
0
 0
+∞
Z
Z
0




−∞

dy,



f (x, y; mT ) dx

f (x, y; mT ) dx




dy,




dy.





(74)

Integrals in the RHS of these equations are readily calculated, since f (x, y; mT )is given by the joint Gaussian probability density in equation (68).
45

1.1.20

Cross–Correlation of 1–Bit Quantized Data:
van Vleck RelationshipIn case of the 1–bit quantization, we have
 x1 = −1,

N = 2, and y1 = −1,

x2 = +1,y2 = +1.

(75)

Therefore, cross–correlation Rx̂ŷ of clipped data x̂[i] and ŷ[i] is given byRx̂ŷ [m] =

2 X
2
X

xi yj P (xi , yj ; m)i=1 j=1=+=−(+1) · (+1) · P (+1, +1; m) + (−1) · (−1) · P (−1, −1; m)(+1) · (−1) · P (+1, −1; m) + (−1) · (+1) · P (−1, +1; m)
P (+1, +1; m) + P (−1, −1; m)
P (+1, −1; m) − P (−1, +1; m).
(76)

Because of symmetric properties of the joint Gaussian probability density
shown in equation (68), the joint probabilities given in equation (74) satisfy
P (+1, +1; m) = P (−1, −1; m),
P (+1, −1; m) = P (−1, +1; m).

(77)

Also, by definition of the probability, sum of probabilities of all possible cases
must be equal to 1, and hence
P (+1, +1; m) + P (−1, −1; m) + P (+1, −1; m) + P (−1, +1; m)= 2 P (+1, +1; m) + 2 P (+1, −1; m) = 1.

(78)

Then, from equations (76), (77), and (78), we obtain
P (+1, +1; m) + P (−1, −1; m)
P (+1, −1; m) − P (−1, +1; m)
2 P (+1, +1; m) − 2 P (+1, −1; m)
4 P (+1, +1; m) − 1.

Rx̂ŷ [m] =
−
=
=

(79)

Substituting the explicit form of the joint Gaussian probability density
in equation (68) into equation (74), we calculte 4 P (+1, +1; m):
4 P (+1, +1; m) = 4

+∞
Z
0

 +∞
Z


0




f (x, y; mT ) dx dy
46

=



 +∞
+∞
Z 
Z
0




0

2
πσx σy

q

x2
y2
xy
+
−
2
r
xy
σx2
σx σy σy21−2
e 2(1 − rxy )2
1 − rxy

!

dx











dy,
(80)

where T is the sampling interval, as before, and we denoted the cross–correlation coefficient rxy [m] = rxy (mT ) as rxy , for simplicity.
Let us introduce new variables ζ and φ, which satisfyx = σx ζ cos φ,
y = σy ζ sin φ.

(81)

Then the above integral is reduced toπ

4 P (+1, +1; m) =

2
π

q

2
1 − rxy
0

π
2

=

Z2

dφ

ζ 2 (1 − rxy sin 2φ)
2 )2(1 − rxy
ζ dζ
e

Z∞ −
0

q

2
1 − rxy2 Z
dφ .
π
1 − rxy sin 2φ

(82)

0

If we further introduce another new variable θ, which satisfiestan φ − rxy 
and, therefore, θ = arctan  q,
2
1 − rxytan φ − rxy,tan θ = q2
1 − rxy

then we have

dθ
=
dφ

q

2
1 − rxy
1
=
.

2 q
2 cos2 φ
1
−
r
sin
2φ
tan φ−rxy
xy
1
−
r
xy
√

1

1

1+



(83)

(84)

2
1−rxy

Note that this has the same form as the integrand of equation (82). The
πlimits of the integration φ = 0 and φ = now correspond to
2


π
r
xy
 and θ = , respectively. Terefore, equationθ = θ0 ≡ − arctan  q
2
2
1 − rxy(82) becomesπ

2
4 P (+1, +1; m) =
π

Z2

θ0





2
2
rxy 
dθ = 1 − θ0 = 1 + arctan  q. (85)2π
π
1 − rxy47

Denoting the cross–correlation coefficient rxy through a sine function:
rxy = sin ξ, and, therefore, q

rxy2
1 − rxy= tan ξ,we obtain
4 P (+1, +1; m) = 1 +2
2
2
arctan(tan ξ) = 1 + ξ = 1 + arcsin(rxy ). (86)πππ

We must specify here a range of arcsin(rxy ), since, in general, arcsine isa multi–value function, while the probability P (+1, +1; m) is certainly not.
In view of the general property of the probability, P (+1, +1, m) must be
confined within a range:
1
0 ≤ P (+1, +1, m) ≤ .
2
Indeed, the upper limit corresponds to a case of the complete correlation(identical data), for which1
P (+1, +1; m) = P (−1, −1; m) = ,2
becauseP (+1, −1; m) = P (−1, +1; m) = 0,
while the lower limit corresponds to a case of the complete anti–correlation(identical data but with different signs), for which
P (+1, +1; m) = P (−1, −1; m) = 0.
Since the cross–correlation coefficient rxy of the original unclipped data mustbe 1 in the complete correlation, and −1 in the complete anti–correlation,
the arcsine function in equation (86) must be confined within a range:−ππ≤ arcsin(rxy ) ≤ .
2
2

(87)

Finally, from equation (79), we obtainRx̂ŷ [m] = 4P (+1, +1; m) − 1 ==2arcsin(rxy [m])π

Rxy (mT ) 2arcsin  q.
πRxx (0) Ryy (0)
48

(88)

Since x̂[i] x̂[i] = 1 and ŷ[i] ŷ[i] = 1 for any i, and therefore Rx̂x̂ [0] = 1 andRŷŷ [0] = 1, for the 1–bit quantized data, their cross–correlation coefficient isequal to their cross–correlation:rx̂ŷ [m] = qRx̂ŷ [m]Rx̂x̂ [0] Rŷŷ [0]= Rx̂ŷ [m].

(89)

Therefore, equation (88) is described also as
rx̂ŷ [m] =2arcsin(rxy [m]).
π(90)

In a particular case of a small cross–correlation coefficient | rxy [m] |  1,which is usually the case in radio interferometry, we have an approximate
linear equation:2(91)
rx̂ŷ [m] = rxy [m].
π2Although we derived equations (88) and (90) assuming that rxy< 1, it isworth to confirm here that the resultant equations are valid in the limitingcases of the complete correlation (rxy = 1 and rx̂ŷ = 1) and the completeanti–correlation (rxy = −1 and rx̂ŷ = −1), too.
Equation (88) or (90) is generally called the “van Vleck relationship”.
This is indeed a surprising result which shows that the cross–correlation
coefficient rxy [m] of the original analog data is almost completely restoredfrom the cross–correlation coefficient rx̂ŷ [m] of the 1–bit quantized data bya simple equation:π
rxy [m] = sinrx̂ŷ [m] .
(92)
2
In the case of the small cross–correlation coefficient | rxy [m] |  1, we haverxy [m] =πrx̂ŷ [m].
2(93)

As we saw before, cross–correlation Rx̂ŷ [m] (= rx̂ŷ [m]) of digital data isreadily obtained from a digital correlator. Therefore, equations (92) and (93)mean that, we can completely derive functional form and, therefore,
spectral shape of the original cross–correlation coefficient from an
output Rx̂ŷ [m] of a digital correlator of the 1–bit quantized data, though theamplitude is reduced by a factor of ∼= 2/π (see, for example, Fifure 21).
It is natural that the cross–correlation Rxy (mT ) itself of the original analog data cannot be directly obtained from the 1–bit quantized data alone
which carries sign information only. Nevertheless, we can restore the cross–correlation Rxy (mT ) of the original analog data, if their dispersions Rxx (0)49Figure 21: Spectra of water maser lines derived from 1–bit quantized data.
Right and left panels show water maser lines of Asymptotic Giant Branch(AGB) stars RT Vir and VY CMa, respectively, obtained at two epochs(top and bottom panels). Thin lines show total–power (i.e. autocorrelation) spectra obtained with Mizusawa 10m antenna. Thick lines show VLBI
cross-power spectra obtained with Mizusawa 10m - Kagoshima 6m baseline.
The cross–power spectra show lower flux densities compared with the total–
power spectra because maser features are slightly extended and hence partially resolved in the VLBI baseline. (Figure courtesy of Imai et al., Astron.
Astrophys, 317, L67-L70, 1997.)

50

and Ryy (0) are known, for example, from suitable system–noise measurements.
Derivation of this relationship owes to the clever use of the probability relations in equation (74) by van Vleck. Note that J.H. van Vleck is ascientist who in 1977 received Nobel Prize in Physics in his major field of
research: the “fundamental theoretical investigations of the electronic structure of magnetic and distorted systems”.
1.1.21

van Vleck Relationship in Autocorrelation

Although we examined the derivation of the van Vleck relationship for cross–correlation, there is no restriction in applying the same logic to autocorrelation of a zero–mean stationary random process x(t). As a matter of fact,
van Vleck originally derived his relationship for the autocorrelation. The
only difference is that we have to use, instead of equation (68), second–order
Gaussian probability density for values of x(t) taken at times t + τ and t:

f (x1 , x2 ; τ ) =1q

2 (τ )
2πσx2 1 − rxx1−2
e 2(1 − rxx (τ ))x21 − 2rxx (τ ) x1 x2 + x22
σx2
,
!

(94)
R
(τ
)
xx
is the correlation
where σx2 = Rxx (0) is the dispersion and rxx (τ ) =Rxx (0)coefficient of x(t).
Exactly the same logic as in the cross–correlation case leads to an equation:!

Rxx (mT )22,rx̂x̂ [m] = Rx̂x̂ [m] = arcsin(rxx [m]) = arcsinππRxx (0)(95)

for the correlation coefficient of the 1–bit quantized data x̂[i], which is thevan Vleck relationship in the autocorrelation.
1.1.22

Spectra of Clipped Data

The derivation of the van Vleck relationship, as we saw above, is quite generaland applicable to any data obeying joint Gaussian probability density. In
particular, the relationship is well valid for continuous–time data, though wehave dealt with only sampled data in our discussion of digitization of analog
data. Thus, we generally haveRx̂ŷ (τ ) =2
arcsin(rxy (τ )),
π
51

(96)

for the cross–correlation, andRx̂x̂ (τ ) =2
arcsin(rxx (τ )),π

and Rŷŷ (τ ) =2
arcsin(ryy (τ )),π(97)

for the autocorrelations, with arbitrary time interval τ .
Therefore, we can calculate cross–power spectrum Sx̂ŷ (ω), and power
spectra Sx̂x̂ (ω) and Sŷŷ (ω) of clipped data in terms of the ordinary Fouriertransformations:Sx̂ŷ (ω) =Sx̂x̂ (ω) =Sŷŷ (ω) =

Z∞

−∞
Z∞

−∞
Z∞

Rx̂ŷ (τ ) e

−iωτ

2
dτ =
π

Rx̂x̂ (τ ) e−iωτ dτ =Rŷŷ (τ ) e−iωτ dτ =

−∞

2
π
2
π

Z∞

−∞
Z∞
−∞
Z∞

arcsin(rxy (τ )) e−iωτ dτ,
arcsin(rxx (τ )) e−iωτ dτ,
arcsin(ryy (τ )) e−iωτ dτ.

(98)

−∞

Let us consider a case when analog data have rectangular spectra. This is
an important case from a practical point of view, since rectangular bandpass
filters are widely adopted in radio astronomy. We saw earlier that correlation
coefficients of analog data take sinc function forms when their spectra arerectangular. Then, how spectra of the data will look like after clipping?
Figure 22 shows a power spectrum Sx̂x̂ (ω) of 1–bit quantized data (solidline in bottom panel) which is derived from analog data with a rectangularpower spectrum Sxx (ω) of bandwidth B (broken line in bottom panel). The
power spectrum of the analog data Sxx (ω) is normalized by the dispersion
Rxx (0) here. In other words, the “power spectrum” in this Figure is theFourier transform of the correlation coefficient rxx (τ ). Therefore, areas under
the spectra (i.e. powers) of the analog and clipped data are equal. Also
shown are (1) correlation coefficient of the analog data rxx (τ ) having a sincfunction form (upper left panel), and (2) correlation coefficient of the clipped2
data rx̂x̂ (τ ) = Rx̂x̂ (τ ) = arcsin(rxy (τ )) (upper right panel).
π
Although the bulk of the spectrum after clipping still remains nearlyband–limited, peak amplitude is somewhat reduced and a low level–skirt appears outside the original band which spreads over a wide range of frequency.
Thus, the original Nyquist interval (1 / (2B)) is no longer strictly optimum
for sampling the clipped data.
On the other hand, Figure 23 shows a cross–power spectrum Sx̂ŷ (ω) of 1–
bit quantized data x̂ and ŷ (solid line in bottom panel) which are derived from52

Correlation coefficient of analog data

Correlation coefficient of clipped data

1

1

0.8

0.8

0.6

0.6

0.4

0.4

0.2

0.2

0

0

-0.2
-10

-0.2
-5

0
2Bτ

5

10

-10

-5

0
2Bτ

5

10

Power spectrum of clipped data
1
0.8
0.6
0.4
0.2
0
-8

-6

-4

-2

0

2

4

6

8

ν/B

Figure 22: Original rectangular power spectrum Sxx (ω) with bandwidth B(broken line in bottom panel) of analog data, and power spectrum Sx̂x̂ (ω)(solid line in bottom panel) of clipped data which is derived from the analog
data by means of the 1–bit quantization. Horizontal axis of the bottom
panel shows frequency ν normalized by the bandwidth B. Upper panels
show correlation coefficient rxx (τ ) of the original analog data having a sincfunction form (left), and correlation coefficient rx̂x̂ (τ ) of the clipped data(right). Horizontal axes of the upper panels show delay τ normalized by the
Nyquist interval (1 / (2 B)) of the original analog data.

53

Cross-correlation coefficient of analog data

Cross-correlation coefficient of clipped data

0.001

0.001

0.00080.00080.00060.00060.00040.00040.00020.000200-0.0002-10-0.0002
-5

0
2Bτ

5

10

-10

-5

0
2Bτ

5

10

Cross-power spectrum of clipped data
0.001
0.00080.00060.00040.0002
0
-8

-6

-4

-2

0

2

4

6

8

ν/B

Figure 23: Original rectangular cross–power spectrum Sxy (ω) with bandwidth B (broken line in bottom panel) of analog data, and cross–power spectrum Sx̂ŷ (ω) (solid line in bottom panel) of clipped data which are derived
from the analog data by means of the 1–bit quantization. Horizontal axis of
the bottom panel shows frequency ν normalized by the bandwidth B. Upper
panels show cross–correlation coefficient rxy (τ ) of the original analog data
having a sinc function form (left), and cross–correlation coefficient rx̂ŷ (τ ) of
the clipped data (right). The maximum cross–correlation coefficient of theanalog data is assumed to be 0.001. Horizontal axes of the upper panels show
delay τ normalized by the Nyquist interval (1 / (2 B)) of the original analog
data.

54

analog data x and y having a real rectangular cross–power spectrum Sxy (ω)
of bandwidth B (broken line in bottom panel). The cross–power spectrum of
the analog data Sq
xy (ω) in this Figure is normalized by the geometric meanof the dispersions Rxx (0) Ryy (0), and, therefore, is the Fourier transform of
the cross–correlation coefficient rxy (τ ). We assume here that the maximum
cross–correlation coefficient of the analog data is small, and, for definiteness,
set it to be 0.001. Upper panels of the Figure show (1) cross–correlation
coefficient of the analog data rxy (τ ) having a sinc function form (left), and (2)2cross–correlation coefficient of the clipped data rx̂ŷ (τ ) = Rx̂ŷ (τ ) ∼= rxy (τ )π(right).
Unlike in the power spectrum case, the spectrum after clipping remains
rectangular, though peak amplitude is reduced by a factor of 2/π ∼= 0.6366
compared with the original analog one, in this cross–power spectrum case.
This is because the cross–correlation coefficient of the clipped data is justproportional to the cross–correlation coefficient of the original analog datawith a proportionality coefficient of 2/π, in the approximation of the small
cross–correlation coefficient | rxy (τ ) | 1, which we assumed here.
1.1.23

Price’s Theorem

Now we proceed to the 2–bit quantization case. Although we can derive
the cross–correlation of the 2–bit quantized data by evaluating probabilitiesof quantization states, as we learned in the 1–bit quantization case, use ofPrice’s theorem (Price, 1958) offers a simpler solution.
The theorem states the following.
Suppose we have zero–mean stationary random processes x(t) and y(t)which satisfy joint Gaussian probability density:f (x, y) ≡ f (x, y; τ ) =2πσx σy1√

1 − r21−2
e 2(1 − r )

x2
xy
y2
− 2r
+
σx2
σx σy σy2 ,where σx2 = Rxx (0) and σy2 = Ryy (0) are dispersions and r ≡ rxy (τ ) =!

Rxy (τ )σx σyis a cross–correlation coefficient of the processes x(t) and y(t).
Then, for an arbitrary function g(x, y) with expectation
hg(x, y)i =Z∞ Z∞g(x, y)f (x, y) dx dy,

−∞ −∞

55

(99)

we have1 ∂ n hg(x, y)i=σxn σyn∂ rn

∂ 2n g(x, y),∂xn ∂y n(100)(Price, 1958).
Auxiliary Formula
In order to prove Price’s theorem, we first derive an auxiliary formula(Papoulis, 1984):∂ 2n f (x, y)
1 ∂ n f (x, y)=,σxn σyn∂r n∂xn ∂y n∂ n f (x0 , y 0 )∂ 2n f (x0 , y 0 )=,∂r n∂x0 n ∂y 0 n

or

(101)

where
x0 =

f (x0 , y 0 ) =

x
,
σx

y0 =
1
√

2π σx σyy,
σyx0 2 − 2 r x 0 y 0 + y 0 22(1 − r 2 )
e
.
2
−

1−r

(102)

We verify this formula using the method of mathematical induction.
• For n = 1, simple differentiations show∂f (x , y )r (1 − r ) + (r x − y ) (r y − x ) −=e∂r
2πσx σy (1 − r 2 )5/2
and
0

0

2

0

0

0

0

x0 2 − 2 r x 0 y 0 + y 0 22(1 − r 2 ),x0 2 − 2 r x 0 y 0 + y 0 2∂ 2 f (x0 , y 0 )r (1 − r 2 ) + (r x0 − y 0 ) (r y 0 − x0 ) −2(1 − r 2 )
=
e
.
0
0
2
5/2
∂x ∂y
2π σx σy (1 − r )∂ 2 f (x0 , y 0 )∂f (x0 , y 0 )=, i.e. the formula is valid for n = 1.
Therefore,∂r∂x0 ∂y 0• Now, if the formula is valid for n = m:
∂ m f (x0 , y 0 )∂ 2m f (x0 , y 0 ),=∂r m∂x0 m ∂y 0 m

56then, for n = m + 1, we have∂ ∂ m f (x0 , y 0 )∂ ∂ 2m f (x0 , y 0 )∂ m+1 f (x0 , y 0 )==∂r m+1∂r∂r m∂r ∂x0 m ∂y 0 m∂ 2m∂ 2(m+1) f (x0 , y 0 )∂ 2m∂f (x0 , y 0 )∂ 2 f (x0 , y 0 ),= 0 m 0m=
= 0m 0m
∂x ∂y
∂r
∂x ∂y
∂x0 ∂y 0∂x0 m+1 ∂y 0 m+1
i.e., the formula is also valid for n = m + 1.
Since the formula is valid for n = 1, this means that the formula isvalid for arbitrary n ≥ 1.
Thus, we confirmed the auxiliary formula.
Proof of Price’s Theorem
We now prove Price’s theorem, using again the method of mathematicalinduction.
• For n = 1, we have1 ∂hg(x, y)i
=
σx σy
∂r
=
=

Z∞ Z∞ "

−∞ −∞
Z∞

−∞

Z∞ Z∞

1 ∂f
g
dx dy =
σx σy ∂r

!

!

2

!

y=+∞

−∞ −∞

∂
∂x

∂f
g
∂y

!

x=+∞

∂f
g
∂y

x=−∞

∂
−
∂y

dy −

Z∞

−∞

Z∞ Z∞

∂2f
dx dy
∂x ∂y

g

−∞ −∞

#

∂g
∂ g
f +
f dx dy
∂x
∂x ∂y
∂g
f
∂x

dx +
y=−∞

Z∞ Z∞

−∞ −∞

∂2g
f dx dy,∂x ∂ywhere we used the auxiliary formula:
1 ∂f (x, y)∂ 2 f (x, y)
=
.
σx σy
∂r
∂x ∂y
2

2

First two terms in the RHS vanish as long as g(x, y) e−(x +y ) converges
to zero at x = ±∞ and y = ±∞. Thus,
1 ∂hg(x, y)i
=
σx σy
∂r

Z∞ Z∞

−∞ −∞

∂ 2 g(x, y)f (x, y) dx dy =∂x ∂yi.e. the theorem is valid for n = 1.
57



∂ 2 g(x, y)
,
∂x ∂y


• Now, if the theorem is valid for n = m:
1 ∂ m hg(x, y)i=σxm σym∂ rm

∂ 2m g(x, y),∂xm ∂y mthen, for n = m + 1, we have1
∂ m+1 hg(x, y)i1 ∂=σxm+1 σym+1∂r m+1σx σy ∂r1 ∂ ∂ 2m g(x, y)==σx σy ∂r ∂xm ∂y m




=

Z∞ Z∞

−∞ −∞

Z∞ Z∞

−∞ −∞

1 ∂ m hg(x, y)iσxm σym∂ rm!

∂ 2m g(x, y) 1 ∂f (x, y)dx dy∂xm ∂y m σx σy∂r∂ 2m g(x, y) ∂ 2 f (x, y)
dx dy.
∂xm ∂y m∂x ∂yIntegrating by parts again, we see that the last term is equal to
Z∞ Z∞

−∞ −∞

∂ 2(m+1) g(x, y)f (x, y) dx dy =∂xm+1 ∂y m+12

∂ 2(m+1) g(x, y),∂xm+1 ∂y m+1
2

as long as g(x, y) e−(x +y ) converges to zero at x = ±∞ and y = ±∞.
Therefore, the theorem is valid for n = m + 1, and, hence, for any
n ≥ 1.
Thus we proved Price’s Theorem.
Price’s Theorem in the 1–Bit Quantaization CaseLet us experience usefulness of Price’s theorem in the 1–bit quantizationcase as an example. If we choose g(x, y) = x̂ ŷ, then we havehg(x, y)i = Rx̂ŷ = hx̂ ŷi,with
x̂(x) =


 +1 :

 −1 :x≥0ŷ(y) =

x < 0,

 +1 :

 −1 :y≥0

y < 0,

in the 1–bit quantization case. We again assume that x and y obey the joint
Gaussian probability density:f (x, y) =2πσx σy1√

1 − r21−2
e 2(1 − r )
58

x2
y2
xy
+
− 2r
σx2
σx σy σy2 .
!

Now, Price’s theorem says
1 ∂Rx̂ŷ
=σx σy ∂r
and we have

∂ 2 (x̂ ŷ)
=
∂x ∂y




∂ x̂ ∂ ŷ
,
∂x ∂y


∂ x̂
∂ ŷ
= 2 δ(x),= 2 δ(y),∂x∂yin the 1–bit quantization case. Therefore∂ x̂ ∂ ŷ∂x ∂yThus we haveand hence



=

Z∞ Z∞

4 δ(x) δ(y) f (x, y) dx dy =−∞ −∞π σ x σy2√

1 − r2.

2∂Rx̂ŷ= √,∂r
π 1 − r22
Rx̂ŷ =
π

Zr
0

√

2
dr 0
arcsin(r),=π
1 − r02where the limits of the integration are chosen to satisfy Rx̂ŷ = 0 at r = 0.
Thus, Price’s theorem allows us to derive the van Vleck relationship in
the really straightforward way.
1.1.24

Cross–Correlation of the 2-Bit Quantized Data

Now we consider, for definiteness, cross–correlation Rx̂ŷ (τ ) of 2–bit quantizeddata x̂(t) and ŷ(t), though logics below are equally applicable to autocorrelation as well. Let clipping criteria for the 2–bit quantization are those given
in Figure 17 and Table 1. We assume a simple case where dispersions of
the original analog data x(t) and y(t) are identical, and denote the common
value as σ, i.e. σx = σy = σ. Then, the joint probability density is given byf (x, y) ≡ f (x, y; τ ) =

2πσ 2

1
√

x2 − 2 r x y + y 2
2σ 2 (1 − r 2 ) ,
e
2
−

1−r

(103)

where r ≡ rxy (τ ) is the cross–correlation coefficient of the analog data.
According to Price’s theorem, the cross–correlation of the clipped data
 2



∂ (x̂ ŷ)∂ x̂ ∂ ŷ
1 ∂Rx̂ŷ=Rx̂ŷ (τ ) satisfies=, where the derivatives
σ 2 ∂r
∂x ∂y
∂x ∂y

59

are now given by∂ x̂
= (n − 1) δ(x + v0 ) + 2 δ(x) + (n − 1) δ(x − v0 ),∂x∂ ŷ= (n − 1) δ(y + v0 ) + 2 δ(y) + (n − 1) δ(y − v0 ).
∂yUsing these derivatives, we can calculate∂ x̂ ∂ ŷ
∂Rx̂ŷ
= σ2
∂r
∂x ∂y




=σ

2

Z∞ Z∞

−∞ −∞

∂ x̂ ∂ ŷ
f (x, y) dx dy.
∂x ∂yAfter simple manipulations, we obtain
(

2

"

2

2

#

)

v
v
v0
1
∂Rx̂ŷ
− 2 0 2− 2 0−
2
σ (1+r) + e σ 2 (1−r)
2σ (1−r ) + 2
= √
+
4(n
−
1)e
(n
−
1)
e
.
∂r
π 1 − r2

Therefore, the solution which satisfies a condition Rx̂ŷ = 0 when r = 0, is
r

2

"

(

2

v
v
1Z
1
− 2 0 0− 2 0 0
2
σ
(1+r
)
σ
(1−r
)
√
(n − 1) e+e
Rx̂ŷ =
π
1 − r020+4(n − 1)e−v2
02σ 2 (1−r 0 2 )

#+ 2 dr 0 .

(104)

In the limiting case where | r | 1, the cross–correlation is given by
"

#

2
2
v0
v0
r
2r
Rx̂ŷ =[(n−1) E +1]2 , (105)
2 (n − 1)2 e− σ2 + 4 (n − 1) e− 2σ2 + 2 =π
πwhere we introduced a notation
2v0

E ≡ e− 2σ2 .

(106)

Note that the cross–correlation of the clipped data Rx̂ŷ (τ ) is proportional tothe analog cross–correlation coefficient r = rxy (τ ), in this case.
1.1.25

Cross–Correlation Coefficient of the 2–Bit Quantized Data

Now we would like to calculate cross–correlation coefficient rx̂ŷ (τ ) of the 2–bit
2
quantized data. For this purpose, we need their dispersions: σx̂x̂= Rx̂x̂ (0)and σŷ2ŷ = Rŷŷ (0), other than the cross–correlation Rx̂ŷ (τ ).
60

We assumed here that analog data x(t) and y(t) have equal dispersions:
Rxx (0) = Ryy (0) = σ 2 , and therefore they obey identical Gaussian probability
densities:
y2
x2
1
1
e− 2 σ2 , and f (y) = √e− 2 σ 2 .
f (x) = √
(107)
σ 2π
σ 2π
Since 4 quantization states in the 2–bit clipped data take values −n, −1,+1, and +n, the joint probabilities of the quantization states of x̂[n] and x̂[n](or ŷ[n] and ŷ[n]), taken at the same time, satisfy
P (−n, −n; 0) + P (−1, −1; 0) + P (+1, +1; 0) + P (+n, +n; 0) = 1, (108)and the dispersion: Rx̂x̂ (0) = hx̂(t) x̂(t)i =4
Xxi xi P (xi , xi ; 0), must be

i=1

Rx̂x̂ (0) = n2 P (−n, −n; 0) + P (−1, −1; 0)+ P (+1, +1; 0) + n2 P (+n, +n; 0).

(109)

Combining equations (108) and (109), and repeating the same thing for ŷ,
we obtainRx̂x̂ (0) = Rŷŷ (0) = Φ + n2 (1 − Φ),(110)
where we introduced a notation:
+v
Z 0

ξ2
1
Φ = P (−1, −1; 0) + P (+1, +1; 0) = √e− 2 σ2 dξ.
σ 2 π −v

(111)0

We used here the explicit form of the probability density in equation (107).
Thus, the cross–correlation coefficient of the 2–bit quantized data rx̂ŷ (τ )is given byrx̂ŷ (τ ) = q=

Rx̂ŷ (τ )Rx̂x̂ (0) Rŷŷ (0)=1

π [Φ + n2 (1 − Φ)]−2v0Rx̂ŷ (τ )
Φ + n2 (1 − Φ)rxy
Z (τ )
0

(

"

2

2

v
v0
1
− 2 0 0−2
σ (1+r ) + e σ 2 (1−r 0 )
√
(n
−
1)
e
1 − r02+4(n − 1)e 2σ2 (1−r0 2 ) + 2



dr 0 ,



#

(112)

as a function of the cross–correlation coefficient of the original analog datar = rxy (τ ). This is an analogue of the van Vleck relationship in the 2–bit
quantization case.
61

In the limiting case of | rxy (τ ) | 1, equations (105) and (110) yieldrx̂ŷ (τ ) =2 [(n − 1) E + 1]2
rxy (τ ).
π [Φ + n2 (1 − Φ)](113)This is an analogue of equation (91) in the 2–bit quantization case.

Correlation coefficient of clipped data

1

1 bit quantization
2 bit quantization

0.8

0.6

0.4

0.2

0

0

0.2

0.4

0.6

0.8

1

Correlation coefficient of analog data

Figure 24: Correlation coefficients of 1–bit (dotted line) and 2–bit (solidline) quantized data with n = 3 and v0 = 0.996 σ as functions of the analog
correlation coefficient.
Exactly the same logic leads to the same functional forms as equations(112) and (113) for a relationship between the correlation coefficient of the
2–bit quantized data rx̂x̂ (τ ) and that of the original analog data rxx (τ ). Thus
the correlation coefficient of the 2–bit quantized data is given byrx̂x̂ (τ ) =1
π [Φ + n2 (1 − Φ)]+4(n − 1)erxx
Z (τ )
0v2− 2 0 022σ (1−r )

(

"

2

2

v
v0
1
− 2 0 0−2
σ (1+r ) + e σ 2 (1−r 0 )
√
e
(n
−
1)
1 − r02+2

dr 0 ,

62

#

(114)

for a general case, andrx̂x̂ (τ ) =2 [(n − 1) E + 1]2
rxx (τ ),π [Φ + n2 (1 − Φ)](115)

for the limiting case of | rxx (τ ) | 1.

Figure 24 shows the relationship between the correlation coefficient, which
could be either (auto–)correlation coefficient or cross–correlation coefficient,of the clipped data and that of the original analog data. The dotted line
shows the van Vleck relationship for the 1–bit quantization case as given byequation (90) or (95). The solid line shows the relationship for the 2–bit
quantization case as given by equation (112) or (114) for a particular set ofparameters n = 3 and v0 = 0.996 σ.
1.1.26

Power and Cross–Power Spectra of 2-Bit Quantized Data

We can now calculate power spectrum Sx̂x̂ (ω) and cross–power spectrumSx̂ŷ (ω)
of the 2–bit quantized data, which are normalized by the dispersion σx̂2 andgeometric mean of the dispersions σx̂ σŷ , respectively, by Fourier transforming
the correlation coefficient rx̂x̂ (τ ) and the cross–correlation coefficient rx̂ŷ (τ )given in equation (112),
Sx̂x̂ (ω) =Sx̂ŷ (ω) =

Z∞

−∞
Z∞

rx̂x̂ (τ ) e

−iωτ

dτ =

rx̂ŷ (τ ) e−iωτ dτ =

Z∞

−∞
Z∞

F (rxx (τ )) e−iωτ dτ,F (rxy (τ )) e−iωτ dτ,

(116)

−∞

−∞

where the function F (r) is given by
(
"
#
Zr
v2
v2
1
1
− 2 0 0− 2 0 0
2
σ
(1+r
)
σ
(1−r
)
√
+e
F (r) =(n − 1) e
02
π (Φ + n2 (1 − Φ)1−r0+4(n − 1)ev2− 2 0 022σ (1−r )

+2





dr 0 .

Let us again consider the case when analog data have rectangular spectra.
Figure 25 shows the power spectrum Sx̂x̂ (ω) of the 2–bit quantized data inthe case with n = 3 and v0 = 0.996 σ (solid line in bottom panel) and that of
the original analog data Sxx (ω) with the rectangular shape of bandwidth B(broken line in bottom panel). Similarly to the 1–bit case, areas under the
spectra (i.e. powers) of the analog and clipped data are equal. Also shownare (1) correlation coefficient of the analog data rxx (τ ) having a sinc function
63

Correlation coefficient of analog data

Correlation coefficient of 2-bit quantized data

1

1

0.8

0.8

0.6

0.6

0.4

0.4

0.2

0.2

0

0

-0.2
-10

-0.2
-5

0
2Bτ

5

10

-10

-5

0
2Bτ

5

10

Power spectrum of 2-bit quantized data
1
0.8
0.6
0.4
0.2
0
-8

-6

-4

-2

0

2

4

6

8

ν/B

Figure 25: Original rectangular power spectrum Sxx (ω) with bandwidth B(broken line in bottom panel) of the analog data, and power spectrum Sx̂x̂ (ω)(solid line in bottom panel) of the 2–bit quantized data with n = 3 and v0 =0.996 σ . Horizontal axis of the bottom panel shows frequency ν normalized
by the bandwidth B. Upper panels show correlation coefficient rxx (τ ) of
the original analog data having a sinc function form (left), and correlation
coefficient rx̂x̂ (τ ) of the 2–bit quantized data (right). Horizontal axes of theupper panels show delay τ normalized by the Nyquist interval (1 / (2 B)) of
the original analog data.

64

Cross-correlation coefficient of analog data

Cross-correlation coefficient of 2-bit quantized data

0.001

0.001

0.00080.00080.00060.00060.00040.00040.00020.000200-0.0002-10-0.0002
-5

0
2Bτ

5

10

-10

-5

0
2Bτ

5

10

Cross-power spectrum of 2-bit quantized data
0.001
0.00080.00060.00040.0002
0
-8

-6

-4

-2

0

2

4

6

8

ν/B

Figure 26: Original rectangular cross–power spectrum Sxy (ω) with bandwidth B (broken line in bottom panel) of the analog data, and cross–power
spectrum Sx̂ŷ (ω) (solid line in bottom panel) of the 2–bit quantized data withn = 3 and v0 = 0.996 σ . Horizontal axis of the bottom panel shows frequency
ν normalized by the bandwidth B. Upper panels show cross–correlation coefficient rxy (τ ) of the original analog data having a sinc function form (left),and cross–correlation coefficient rx̂ŷ (τ ) of the 2–bit quatized data (right).
The maximum cross–correlation coefficient of the analog data is assumed tobe 0.001. Horizontal axes of the upper panels show delay τ which is normalized by the Nyquist interval (1 / (2 B)) of the original analog data.

65

form (upper left panel), and (2) correlation coefficient of the 2–bit quantizeddata rx̂x̂ (τ ) (upper right panel).
The spectrum after clipping again shows somewhat reduced peak amplitude and a low–level skirt extending over a wide range of frequency. Thus,
the original Nyquist interval (1 / (2B)) is not strictly optimum for sampling
the clipped data again, though to a smaller extent compared with the 1–bitcase.
On the other hand, Figure 26 shows a cross–power spectrum Sx̂ŷ (ω) of the
2–bit quantized data x̂ and ŷ in the case with n = 3 and v0 = 0.996 σ (solidline in bottom panel) and that of the original analog data x and y, having areal rectangular cross–power spectrum Sxy (ω) of bandwidth B (broken line inbottom panel). We assumed here again that the maximum cross–correlation
coefficient of the analog data is as small as 0.001, for definiteness. Upper
panels of the Figure show (1) cross–correlation coefficient of the analog datarxy (τ ) having a sinc function form (left), and (2) cross–correlation coefficientof the clipped data rx̂ŷ (τ ) (right). Note that the spectrum after clipping
remains again rectangular, because of the approximate linearity between theanalog and clipped cross–correlation coefficients.
1.1.27

Dispersion of Digital Correlator Output

Let us now examine the signal–to–noise ratio of an output of a digital correlator. How does the S/N ratio differ from the one expected from the analog
or “semi–analog” correlator, which we saw earlier? For this purpose, we firstconsider dispersion of the digital correlator output.
The output R of a digital correlator, which averages products of N samples of clipped data x̂[i] and ŷ[i] from two antennas in a radio interferometer,
is
N
1 X
R=
x̂[i] ŷ[i].
(117)
N i=1

We again assume here that the delay tracking and the fringe stopping arecompletely performed before the multiplication, so that signal parts in thetwo data streams x̂[i] and ŷ[i] are perfectly aligned.
Under the assumption of the stationary random processes, the expecta2tion hRi and the dispersion σRof the output R are given by
hRi = hx̂[i] ŷ[i]i = Rx̂ŷ [0],2σR
= hR2 i − hRi2 =

(118)

N X
N
1 X
2
hx̂[i] ŷ[i] x̂[j] ŷ[j]i − Rx̂ŷ[0]. (119)N 2 i=1 j=166

Note that usually N is a huge number. For example, if we integrate data
sampled with a rate of 4 Msps (mega sample per second) for 1 second, thenN ≈ 4, 000, 000.
Let us assume that the original analog data x(t) and y(t) have rectangularbaseband spectra with bandwith B, and the data are sampled with Nyquist
rate 2B. In such a case, different sample pairs are independent in the analog
data, as we discussed earlier, i.e. Rxy [m], Rxx [m], and Ryy [m] are all zero,
if m 6= 0. The same statement is valid for the 1–bit quantized and 2–bit
quantized data, since correlation coefficients, and correlations, too, of theseclipped data are equal to zero when those of the original analog data arezero, as we can easily see in equations (90), (95), and (112). Thus, Rx̂ŷ [m],Rx̂x̂ [m], and Rŷŷ [m] are also all zero, if m 6= 0. In addition, we assume thatthe cross–correlation Rx̂ŷ [0] is much smaller than the autocorrelations Rx̂x̂ [0]and Rŷŷ [0], as usually so in radio interferometers.
Then in the double sum of equation (119):N X
N
X

i=1 j=1hx̂[i] ŷ[i] x̂[j] ŷ[j]i =

N
X
i=1

hx̂[i] ŷ[i] x̂[i] ŷ[i]i +N X
NX
i=1 j6=ihx̂[i] ŷ[i] x̂[j] ŷ[j]i,(120)dominating terms will be
N X
NX

hx̂[i] ŷ[i] x̂[j] ŷ[j]i ∼=i=1 j=1NX
i=1hx̂[i] x̂[i]i hŷ[i] ŷ[i]i +N X
NX
i=1 j6=ihx̂[i] ŷ[i]i hx̂[j] ŷ[j]i2
2= N Rx̂x̂ [0] Rŷŷ [0] + N (N − 1) Rx̂ŷ
[0] ∼
[0],
= N Rx̂x̂ [0] Rŷŷ [0] + N 2 Rx̂ŷ(121)where we neglected 1 compared with the large number N .
Now equations (119) and (121) yeild an approximate formula for thedispersion of the digital correlator output:
2 ∼
σR
=

1.1.281
1
2
2
Rx̂x̂ [0] Rŷŷ [0] + Rx̂ŷ[0] − Rx̂ŷ[0] =Rx̂x̂ [0] Rŷŷ [0].
N
N

(122)

S/N Ratio of Digital Correlator Output

We are now ready to calculate the signal hRi to noise σR ratio SN R of thedigital correlator output, using equation (122). We obtain
SN R =√√ rx̂ŷ [0]hRi √Rx̂ŷ [0]= N rx̂ŷ [0] = N
= N q
ρ,
σR
rxy [0]Rx̂x̂ [0] Rŷŷ [0](123)

where ρ ≡ rxy [0] is the maximum cross–correlation coefficient of the original
analog data under our assumption of the perfect delay tracking and fringe67

stopping. As we saw earlier, the maximum cross–correlation coefficient ρ fora continuum spectrum source is approximately given byρ =s

TA 1 TA 2,
TS 1 TS 2

where TA1 , TA2 and TS1 , TS2 are antenna temperatures and system noise
temperatures, respectively, at antennas 1 and 2, which are both assumed
constant throuout the frequency band B. On the other hand, if the integration time is τa , the number of samples N with the Nyquist interval 1 / (2 B)
is equal to
N = 2 B τa .
Therefore, equation (123) is reduced torx̂ŷ [0]SN R =
rxy [0]s

TA 1 TA 2 q
rx̂ŷ [0]SN R analog ,2 B τa =TS 1 TS 2
rxy [0]where
SN R analog ≡

s

(124)

TA 1 TA 2 q
2 B τa ,TS 1 TS 2

is the signal–to–noise ratio of an analog correlator output for the same continuum spectrum source received with the same antenna–receiver systems,
as we saw in equation (65).
Thus the so–called “coherence factor” ηc , which is defined as
ηc ≡

SN R
,
SN Ranalog(125)rx̂ŷ [0],rxy [0](126)

is given byηc =in this case.
For the case of the small cross–correlation coefficient | rxy [0] |  1, which
is usually the case in radio interferometry, rx̂ŷ [0] of the clipped data is proportional to rxy [0] of the original analog data. Then, in view of equations(91) and (113), we have2
ηc = ∼(127)= 0.64,π
for the 1–bit quantization case, andηc =2 [(n − 1) E + 1]2,
π [Φ + n2 (1 − Φ)]
68

(128)

for the 2–bit quantization case.
With the coherence factor ηc , equation (124) is reduced to
SN R = ηcsTA 1 TA 2 q
2 B τa ,TS 1 TS 2(129)

which is an important formula for estimating sensitivity of a radio interferometer observing a continuum spectrum source. Note that the coherence
factor as given by equation (126) still takes into account only the loss due
to effects of clipping. Usually, the factor is further reduced in view of losses
occuring through digital processing in hardware correlators.
1.1.29Optimum Parameters v0 and n in the 2–Bit Quantization

We mentioned earlier that the two parameters in the 2–bit quantization, thethreshold v0 and the higher–level value n, are chosen so that the signal–to–noise ratio of the digital correlator output is maximized. This condition is
Coherence factor of 2-bit quantized data

0.885
0.88
0.875
0.87
0.865
0.86

2.5

3
n

3.5

4

0.95
0.9
0.85
4.50.8

1

1.15
1.1
1.05
v0 / σFigure 27: Coherence factor of the 2–bit quantized data as a function of
normalized threshold v0 /σ and higher–level value n in the limiting case of
the small cross–correlation coefficient.
69

fulfilled when the coherence factor ηc takes the maximum. In the limiting caseof the small cross–correlation coeficient | rxy [0] |  1, which is common inradio interferometry, we can calculate the coherence factor ηc using equation(128) with explicit forms for the functions E and Φ given in equations (106)and (111).
Figure 27 shows the coherece factor as a function of the two parameters, v0
normalized by the analog dispersion σ, and n, in the small cross–correlation
coefficient case.
As we see from the Figure, the maximum value of the coherence factor
ηc = 0.883 is obtained when n = 3.34 and v0 = 0.982 σ. However, froma practical point of view, designing of digital circuitry for 2–bit hardware
correlators can be more easily implemented when n is an integer. Therefore,n = 3, v0 = 0.996 σ with ηc = 0.881, and n = 4, v0 = 0.942 σ with ηc = 0.880,
are often used in existing VLBI 2–bit quantization systems. This is why weadopted n = 3, v0 = 0.996 σ in Figures 24, 25, and 26.
Accordingly, in the limiting case of the small cross–correlation coefficient| rxy [0] |  1, the coherence factor, which takes into acount the effect ofclipping only, is given by
ηc ∼
(130)
= 0.64,for the 1–bit quantization case, andηc ∼= 0.88,(131)

for the 2–bit quantization case.
1.1.30Effect of Oversampling in S/N Ratio of Clipped Data

As we saw earlier, power spectra of clipped data show low–level but broad
skirt beyond band edges at ν = ±B of the original rectangular spectra ofcorresponding analog data (Figures 22 and 25). Therefore, the Nyquist sampling for the original analog data with 2B rate is no longer optimum for theclipped data.
Sampling with a rate higher than the Nyquist rate for the analog data 2B,
which is a little incorrectly called the “oversampling”, improves the signal–
to–noise ratio of the clipped data. In this case, different sample points are nolonger independent any more, and we have to take into account contributionsof autocorrelations between different sample points, i.e. Rxx [m] and Ryy [m]with m 6= 0, when we calculate the signal–to–noise ratio.
Let us consider that we sample our data with a rate which is β times as
fast as the Nyquist rate of the original band–limited analog data 2B. Then,in the calculation of the dispersion of the digital correlator output R shown70

in equations (119) – (122), we must leave autocorrelation terms with non–zero arguments. Specifically, the second term in the RHS of equation (120)now must be
N X
NX
i=1 j6=ihx̂[i] ŷ[i] x̂[j] ŷ[j]i ∼
=

N X
N
X
i=1 j6=i[hx̂[i] ŷ[i]i hx̂[j] ŷ[i]i + hx̂[i] x̂[j]i hŷ[i] ŷ[i]i],where the products of the autocorrelations in the second term in the RHS
could be well larger than the products of the cross–correlations in the firstterm in this oversampling case. Thus the dispersion now becomes
2
σR
=



N
X



N X
N
X

1 
2
hx̂[i] ŷ[i] x̂[i] ŷ[i]i +hx̂[i] ŷ[i] x̂[j] ŷ[j]i − Rx̂ŷ
[0]
2
N i=1
i=1 j6=i

N X
N
1
1 X
∼
Rx̂x̂ [0] Rŷŷ [0] + 2
Rx̂x̂ [i − j] Rŷŷ [i − j],=N
N i=1 j6=i(132)

where the second term in the RHS stands for the largest contribution of
the autocorrelations between different i–th and j–th sample points. Since
number of combinations of i and j having the same diffference k = i − j
is N − | k |, and autocorrelations are likely to be large enough for small| k |  N only, we have
2 ∼
σR
=

−1
2 NX
1
Rx̂x̂ [0] Rŷŷ [0] + 2(N − k)Rx̂x̂ [k] Rŷŷ [k]
N
N k=1
∞
X

1
∼
Rx̂x̂ [0] Rŷŷ [0] (1 + 2
rx̂x̂ [k] rŷŷ [k]).
=
N
k=1

(133)

Thus, for a continuum spectrum source, the signal–to–noise ratio of thedigital correlator output of N = 2 β B τa oversampled data is approximately
given bys√hRirx̂ŷ [0] TA1 TA2
2 β B τav
SN R =
=
.
(134)
u
∞
σR
rxy [0] TS1 TS2 u
X
t1 + 2rx̂x̂ [k] rŷŷ [k]k=1

Therefore, the coherence factor in the case of the oversampling is given by√
SN Rrx̂ŷ [0]
β
v
ηc =
.
(135)
=
∞
SN Ranalog
rxy [0] uXut1 + 2rx̂x̂ [k] rŷŷ [k]k=171We can estimate this factor, using equations (95) and (114) which describe
correlation coefficients of the clipped data rx̂x̂ [k] and rŷŷ [k] as functions ofanalog correlation coefficients rxx [k] and ryy [k], which we denote as rx̂x̂ [k] =fcl (rxx [k]) and rŷŷ [k] = fcl (ryy [k]). Here fcl is the clipping function which is
given in equation (95) and equation (114) for the 1–bit and 2–bit quatization
cases, respectively.
When the original analog data have rectangular baseband spectra withbandwidth B, as we have assumed in the present discussion, the analog
correlation coefficients rxx (τ ) and ryy (τ ) have sinc function forms:
sin(2π B τ ),2π B τ
as given in equation (46). Therefore, correlation coefficients of their timesamples x[k] and y[k] with the sampling interval Ts = 1/(2 B β) are given byrxx (τ ) = ryy (τ ) =rxx [k] = ryy [k] = rxx (k Ts ) = ryy (k Ts ) =

sin



πk
β
πk
β



.

(136)

Coherence Factor vs Oversampling Ratio
1
0.8

ηc

0.6
0.4
0.2
0

analogtwo-bit
one-bit
0

2

4

6

8

10
β

12

14

16

18

20

Figure 28: Coherence factor ηc for the analog (+), 2–bit quantized (×), and1–bit quantized (∗) data as functions of the oversampling ratio β.
Thus, the signal–to–noise ratio and the coherence factor of the oversampled clipped data is given by
SN R = ηcsTA 1 TA 2 q
2 B τa ,TS 1 TS 2
72

(137)

√
β
rx̂ŷ [0]
v
ηc =
 .

u
rxy [0] u
πk
∞
sin
X
u
β 
2 
t1 + 2
fcl
πk
k=1

(138)

β

Figure 28 shows the coherence factor ηc for the analog, 2–bit quantized,
and 1–bit quantized data as functions of the oversampling ratio β, calculated
by means of equation (138). The coherence factor of the clipped data improves gradually with increasing β, approaching to some constant values. On
the other hand, the coherence factor of the analog data is always 1, irrespective of the oversampling ratio β, in accordance with the sampling theorem.
In a particular case of the “double–speed sampling” β = 2, we have ηc =0.744 for the 1–bit quantized data, and ηc = 0.930 for the 2–bit quantizeddata.
1.1.31

Coherence Factor and Sensitivity with Given Bit Rate

Coherence factor ηc with different clipping and oversampling is summarizedin Table 2. It is obvious that the coherence factor increases as we increase
Number of

Number ofcoherence factor ηcbits Nb

quantization levels

β=1

β=2

1

2

0.64

0.74

2

4

0.88

0.93

Table 2: Relationship of coherence factor ηc with number of bits Nb andoversampling ratio β.
the number of bits Nb and the oversampling ratio β.
However, actual observations are usually limited by maximum bit–rateνb of data streams allowed by correlators or recorders. If the maximum bit–rate νb is fixed, the maximum permissible sampling rate, which is β times aslarge as the Nyquist rate 2B with the analog bandwidth B, is νb /Nb = 2 β B.
Therefore, the maximum allowable bandwidth is limited by the bit–rate as
B = νb /(2 β Nb ). Hence, in view of equation (137), the maximum signal–
to–noise ratio, i.e. the sensitivity of the interferometric observation, for acontinuum spectrum source, is proportional to√ηc
SN R ∝ ηc B ∝ √
.
(139)
βNb
73

√
η c / β Nb

Number of

Number ofbits Nb

quantization levels

β=1

β=2

1

2

0.64

0.52

2

4

0.62

0.47

√
Table 3: Factor ηc / β Nb , which determines the signal–to–noise ratio, estimated for clipped and oversampled data with a fixed bit–rate.
√
Table 3 shows the factor ηc / β Nb , which determines the signal–to–noise
ratio for the continuum spectrum source. When the bit–rate is given, thesensitivity turns out to be largest in the simplest case of the 1–bit quantization (Nb = 1) with the Nyquist sampling (β = 1)! This is the reason why
many VLBI observations still use the 1–bit quantization scheme with the
Nyquist sampling.
However, if we observe a line spectrum source, such as an astronomical
maser source, the frequency range containing the signal from the radio sourceis confined within a limited spectral profile of the source where we often see
many narrow lines. Thus, in order to estimate the signal–to–noise ratio ofa spectral line, we must replace in equation (137) the bandwidth B of ourreceiving system by the width of the spectral line. The line width is intrincic
to the radio source and constant irrespective of system–dependent parameters, such as the number of bits Nb , or the oversampling ratio β. Therefore,
the larger the coherence factor ηc , the higher is the sensitivity in this case.
Consequently, modern VLBI systems tend to adopt the 2–bit quantization
scheme for better performance in line–spectrum–source observations. Notealso that the 2–bit quantization scheme (Nb = 2) with the Nyquist sampling(β = 1) offers almost the same sensitivity as the 1–bit quantization schemewith the Nyquist sampling, i.e. 0.62 against 0.64, for a continuum source,
as evident from Table 3. Figure 29 shows an example of the maser source
spectra derived from 2–bit quantized data.

74

Figure 29: Spectra of the water maser source Cep A in an active region ofstar formation derived from 2–bit quantized data. Upper panel shows total–
power spectrum with Mizusawa 20m antenna, and lower panel shows phaseand amplitude of cross–power spectrum with Mizusawa–Ogasawara baselineof the VERA. Profiles of the two spectra are fairly different which indicate
effects of partial resolution of maser features in the VLBI baseline. Notethat amplitude scales are different in the two panels. (T. Hirota, private
communication in 2005.)
75

1.2

Frequency Standard

1.2.1

VLBI Requires “Absolute” Frequency Stability
RF

ν

IF

ν - ν0analogdigitalConnected Interferometer
instrumental delay
Multiplier

Integrator

correlator

Tape recorderMultiplier

VLBI

Integrator
correlator

Figure 30: Connected–Element Interferometer vs. VLBI system.
Development of the highly stable frequency standard was indispensable
for realization of VLBI.
In the case of the connected–element interferometer, a common frequency
standard is used to generate local oscillator (LO) reference signals for frequency conversion in all antennas (see upper panel of Figure 30). Therefore,
any phase noise, due to instability in the frequency standard, is common in
the signals from all antennas. Such common phase noise from any given antenna is always compensated by the same phase noise from another antenna,
in the correlation processing.
In fact, let us assume that a signal from a radio source, which we assume
to be point–like for simplicity, received by antenna A is:vA (t) ∝ sin[ωt − φ(t)],and the same signal received by antenna B is:vB (T ) ∝ sin[ω(t − τg ) − φ(t)],76where ω is a frequency of the radio signal, τg is the geometric delay between
the two antennas, and φ(t) is the common phase noise due to the frequency
standard.
After multiplication of these signals, we havevA (t) vB (t) ∝1{cos(ωτg ) − cos[2 ωt − ωτg − 2 φ(t)]},2

and after integration (averaging), we havehvA (t) vB (t)i ∝1cos(ωτg ),
2

(140)

since rapidly oscillating second term with frequency 2 ω is averaged out.
Thus, no effect of the phase noise remains in the correlator output! Consequently, correlation results in connected–element interferometers
are almost unaffected by the instability of the frequency standard.
Then, we readily obtain the almost pure fringe pattern ∝ cos(ωτg ) in thecrorrelation results as shown in equation (140), as far as we perform sufficiently long integration to suppress the thermal noise occuring in receiving
systems and in the environment, and achieve high enough signal–to–noise ratio. This is why very high stability of frequency standards is not necessarily
required in connected–element interferometers.
How about VLBI, then?
Each antenna in VLBI uses its own independent frequency standard togenerate the LO reference signal (see lower panel of Figure 30).
Then, a signal from a point–like radio source received by antenna A isgiven byvA (t) ∝ sin[ωt − φA (t)],and the same signal received by antenna B is:vB (t) ∝ sin[ω(t − τg ) − φB (t)],where φA (t) and φB (t) are phase noises in independent frequency standards.
After multiplication, we have
1vA (t) vB (t) ∝ {cos[ωτg − φA (t) + φB (t)] − cos[2 ωt − ωτg − φA (t) − φB (t)]},2
and after integration (averaging),
hvA (t) vB (t)i ∝1cos[ωτg − φA (t) + φB (t)].
2
77

(141)

Therefore, no compensation of the phase noises is expected in the correlator
output of VLBI. In other words, correlation results in VLBI are always
directly affected by the instability of the frequency standards!
The noise in the fringe phase (the argument of the sinusoidal fringe pattern in equation (141)) gives rise to two difficulties in VLBI.
First, the fringe phase, which is the important observable in radio interferometry, is contaminated by the phase noise.
Second, the phase noise severely limits the sensitivity of VLBI. Indeed, it
becomes impossible for us to completely stop the oscillation of the sinusoidal
fringe pattern if the phase noise varies in time, even when we ideally compensate the geometric delay τg by applying accurate enough delay tracking andfringe stopping. Then, if we further integrate (time–average) the correlator
output, hoping to get higher signal–to–noise ratio, the phase–noise–induced
oscillation of the fringe pattern results in smaller amplitude of the averaged
signal. The thermal noise contribution in the correlator output must be
surely suppressed by the integration, but the averaged signal itself could be
reduced even more rapidly due to the oscillation. Depending on the ratio
between the integration time and the timescale of the oscillation, the integration may not improve the signal–to–noise ratio at all, but even degrade
it. Such an effect is called the “coherence loss”. Therefore, the integration
time must be short enough not to cause large coherence loss, but this implies
that the signal–to–noise ratio must be limited by the short integration time.
This is why we need “absolute” stability of frequency standards in VLBI,
in order to ease these difficulties.
Of cource, not only the noise due to the frequency standards, but anyother phase noise, due for example to the propagation delay through the
turbulent atmosphere, causes similar difficulties, as we will see later.
Now we would like to discuss the way to quantitatively describe the frequency stability, and to estimate effects of the incomplete stability.
1.2.2

How to Describe Frequency Stability?

Let us consider a reference signal v(t) which is expected to have a form in
an ideal case:v(t) = v0 cos(ω0 t),(142)
with a nominal frequency ω0 . In actuality, however, any real reference signal
has a phase noise φ(t), and therefore has a form:v(t) = v0 cos(ω0 t + φ(t)).

78

(143)

In this actual case, instantaneous frequency ωa (t) will beωa (t) = ω0 +

dφ(t)
.
dt

(144)

We then introduce a concept of the “fractional frequency deviation(FFD)” y(t), which is defined byy(t) =δω(t)ωa (t) − ω0
1 dφ(t)
=
=
,
ω0
ω0
ω0 dt

(145)

as a measure of the frequency stability.
Let us assume that φ(t) and y(t) are stationary random processes. Then,
their autocorrelations are functions of the time difference τ :Rφφ (τ ) = hφ(t + τ ) φ(t)i, and Ryy (τ ) = hy(t + τ ) y(t)i,(146)

and their power spectra are given by Fourier transform relations:Sφφ (ω) =Syy (ω) =

Z∞

−∞
Z∞

Rφφ (τ ) e

−iωτ

dτ,

Ryy (τ ) e−iωτ dτ,1
Rφφ (τ ) =2π
Ryy (τ ) =

−∞

1
2π

Z∞

−∞
Z∞

Sφφ (ω) eiωτ dω,Syy (ω) eiωτ dω. (147)−∞

In view of equation (145), the autocorrelations of y(t) and φ(t) are mutually related by1
Ryy (t, t ) = hy(t) y(t )i = 2
ω0
0

0



dφ(t) dφ(t0 )
dt
dt0



=

1 ∂2
Rφφ (t, t0 ). (148)ω02 ∂t ∂t0

In our case of the stationary random processes, this is reduced toRyy (τ ) = −1 d2
Rφφ (τ ),ω02 dτ 2(149)

where τ = t − t0 .
From equations (147) and (149), we have1 d2 1
Ryy (τ ) = − 2 2ω0 dτ 2πZ∞

Sφφ (ω) e

iωτ

−∞

1
dω =
2π

Z∞

−∞

ω2
S (ω) eiωτ dω,2 φφ
ω0and, therefore, the power spectra of y(t) and φ(t) are related to each otherby a relation:
ω2
Syy (ω) = 2 Sφφ (ω).
(150)
ω0
79

Let us make a comment here about conventions which have been traditionally used for describing power spectra in frequency stability discussions.
Since y(t) and φ(t) are real functions of time, the autocorrelation Rφφ (τ )and Ryy (τ ), and power spectra Sφφ (ω), Syy (ω) are all real and even functionsof τ and ω. Using this property, and using frequency ν instead of angularfrequency ω = 2π ν, we can describe the power spectra in the “single–sided
forms” Sφφ (ν) and Syy (ν), which have been widely used in the frequency
stability discussions:Rφφ (τ ) =

0

Z∞

Z∞

Ryy (τ ) =Z∞

Sφφ (ν) = 4

Z∞

Syy (ν) = 4

Rφφ (τ ) cos(2π ν τ ) dτ,Sφφ (ν) cos(2π ν τ ) dν,0

Ryy (τ ) cos(2π ν τ ) dτ,Syy (ν) cos(2π ν τ ) dν.

0

0

(151)
These single–sided spectra Sφφ (ν) and Syy (ν) are related to our double–sided
power spectra Sφφ (ω) and Syy (ω) by:
Sφφ (ν) = 2 Sφφ (2πν), and Syy (ν) = 2 Syy (2πν),(152)

for the positive frequency range ν ≥ 0. Therefore, their mutual relationshiphas been often given by:
Syy (ν) =ν2
Sφφ (ν),ν02

instead of equation (150), where ν0 =

(153)

ω0
.
2π

However, we will continue to use our double–sided power spectrum forms
in equations (147), following our previous discussions.
1.2.3

Types of Phase and Frequency Noises

Measurements have shown that noises in the frequency stability are classified into following types according to the power–law index α of the power
spectrum of the FFD Syy (ω) = Hα ω α with roughly constant coefficients Hα .
Each noise type, or power–law component of the FFD spectrum, has its own
characteristic name listed in Table 4.
An oscillator in a frequency standard often shows a combination of all
power–law components of Table 4 in various frequency ranges, as shown
schematically in Figure 31.
80

αName of Noise TypeSyy (ω)Sφφ (ω)2

White phaseH2 ω 2

ω02 H2

1

Flicker phaseH1 ω 1

ω02 H1 ω −10White frequency

H0

ω02 H0 ω −2

-1

Flicker frequency

H−1 ω −1

ω02 H−1 ω −3

-2

Random walk of frequency

H−2 ω −2

ω02 H−2 ω −4

Table 4: Power–law types of phase and frequency noise spectra.

α=−2

α=2

α=−1

α=1
α=0

Figure 31: A schematic spectrum of the FFD, y(t), of an oscillator showing
all power–law components in Table 4 in various frequency ranges. Horizontal
and vertical axes show the angular frequency ω and the power spectrum
Syy (ω), respectively, in log–scales.

811.2.4

Time Domain Measurements

It has traditionally been easier to make measurements in the time domain,than in the frequency domain. Hence, frequency stability characteristics areusually given in terms of the time domain measurements.
Suppose that we have N phase values: φ(t1 ), φ(t2 ), · · · , φ(tk ), · · · , φ(tN ),measured at equally spaced time points: t1 , t2 , · · · , tk , · · · , tN , with timeinterval T , where T = tk+1 − tk for any k.
Then, using these phase values measured in the time domain, we form a
discrete time series of the “mean fractional frequency deviation” ȳ[k]which is defined by
tZk+1

1
φ(tk+1 ) − φ(tk )=ȳ[k] =
ω0 T
ω0 T

tk

dφ(t0 ) 0
1
dt =
0
dt
T

tZk+1

y(t0 )dt0 ,

(154)

tk

and can be described also as
1ȳ[k] =
T

t+ T2

Z

y(t0 )dt0 ,

(155)

t− T2

where t = tk + T /2. This is a running mean of the FFD y(t) at the pointt = tk + T /2 over the time interval T . Thus, ȳ[k] for any k is given by alinear system of y(t):ȳ[k] = y(t) ∗ a(t) =Z∞−∞y(t − t0 ) a(t0 ) dt0 ,(156)

with an impulse response a(t):a(t) =





1
T


 0

if − T2 < t ≤ T2 ,otherwise,(157)

at t = tk + T /2, similarly to what we saw in Chapter 3.
1.2.5“True Variance” and “Allan Variance” of Fractional Frequency Deviation

Using the mean fractional frequency deviation ȳ[k], which is obtained byaveraging the FFD y(t) for the time interval T , as shown in equation (155),82

we introduce the “true variance of the fractional frequency deviation(TVAR)” as a function of the time interval T :I 2 (T ) = hȳ 2 [k]i,(158)

where, we assumed zero–mean of the mean FFD ȳ[k]:hȳ[k]i = 0.
In view of the ergodicity, we can estimate the TVAR I 2 (T ) by means of atime–average Σt of N − 1 values of the squared mean FFD ȳ 2 [k]:−1
1 NXȳ 2 [k].
Σt =N − 1 k=1

(159)

Generally speaking, there seems no difficulty in calculating the time–
average Σt , and, therefore, this estimated TVAR appears a good measureof the frequency stability. In actuality, however, Σt , as given in equation(159), diverges for some important types of frequency noises, and therefore
cannot be used for characterizing the frequency stability as a whole. The
assumption of the stationary random process, or the ergodicity, is not likely
to be strictly fulfilled in these diverging cases.
In order to overcome this difficulty, David W. Allan proposed to take adifference of successive two samples:∆ȳ[k] = ȳ[k + 1] − ȳ[k],(160)

and calculate an average of the squared difference (∆ȳ[k])2 divided by two:N
−2X
1(∆ȳ[k])2 ,Σa =2 (N − 2) k=1(161)(Allan, 1966). In view of equation (154), the two sample difference ∆ȳ[k] ofthe mean FFD is easily derived from the measured phase noise φ[k] by∆ȳ[k] =φ(tk+2 ) − 2 φ(tk+1 ) + φ(tk )
.
ω0 T

(162)

The mean square sum Σa gives an estimation of the “two–sample variance” or “Allan variance (AVAR)”:σy2 (T ) =h(∆ȳ[k])2 i
,
2
83

(163)

y[k]

τ
∆y[k]=y[k+1]-y[k]τFigure 32: Mean square sum of the two–sample difference of the mean FFD(lower panel) has a better chance to converge than that of the mean FFD
itself (upper panel), even when the mean FFD shows a diverging long–time–scale behaviour.
as a function of the time interval T .
Estimation of the AVAR with Σa has a better chance to converge thanthat of the TVAR with Σt , as schematically illustrated in Figure 32. Therefore, the estimated AVAR is widely accepted as a good measure of the frequency stability.
In view of equation (154), we can describe the two–sample difference ofthe mean FFD ∆ȳ[k] through time integrations of the FFD y(t):∆ȳ[k] = ȳ[k + 1] − ȳ[k] == tZk+2
1 

T 


tk+1


 t+T
Z
1 

T 


t

y(t0 ) dt0 −00y(t ) dt −

tZk+1
tk

Zt





y(t0 ) dt0 y(t0 ) dt0 ,

t−T

(164)




where t = tk+1 . This is nothing but a linear system of the FFD y(t):∆ȳ[k] = y(t) ∗ b(t) =Z∞−∞y(t − t0 ) b(t0 ) dt0 ,(165)

with an impulse response:







1
T
b(t) = − T1







0

if − T < t ≤ 0,
if 0 < t ≤ T ,otherwise,
84

(166)

at time t = tk+1 .
1.2.6

True Variance and Allan Variance through Power Spectrumof Fractional Frequency Deviation

Equations (156) and (165) describe the discrete–time series ȳ[k] (the meamFFD) and ∆ȳ[k] (the two–sample difference of the mean FFD) through thelinear systems with impulse responses a(t) and b(t) at times t = tk + T /2
and t = tk+1 , respectively. We can formally extend these linear systems toyield continuous–time functions ȳ(t) and ∆ȳ(t):ȳ(t) = y(t) ∗ a(t) =∆ȳ(t) = y(t) ∗ b(t) =

Z∞

−∞
Z∞

−∞

y(t − t0 ) a(t0 ) dt0 ,y(t − t0 ) b(t0 ) dt0 ,(167)

at arbitrary time t, and introduce their autocorrelations
Rȳȳ (τ ) = hȳ(t + τ ) ȳ(t)i,R∆ȳ∆ȳ (τ ) = h∆ȳ(t + τ ) ∆ȳ(t)i,(168)

and power spectra Sȳȳ (ω) and S∆ȳ∆ȳ (ω). Of course, ȳ(t) takes the particularvalue ȳ[k] at t = tk + T /2, and ∆ȳ(t) takes the particular value of ∆ȳ[k] at
t = tk+1 , i.e.
T,2∆ȳ(t) = ∆ȳ[k] at t = tk+1 .
ȳ(t) = ȳ[k]at t = tk +(169)

Then, we can describe the TVAR, I 2 (T ) , through the autocorrelationRȳȳ (τ ), and then the power spectrum Sȳȳ (ω):∞

1 Z
Sȳȳ (ω) dω.
I (T ) = h(ȳ[k]) i = Rȳȳ (0) =
2π
2

2

(170)

−∞

Similarly, for the AVAR, σy2 (T ) , we haveR∆ȳ∆ȳ (0)h(∆ȳ[k])2 i
σy2 (T ) =
=
=
2

2

1
4π

Z∞

S∆ȳ∆ȳ (ω) dω.

(171)

−∞

Since the functions ȳ(t) and ∆ȳ(t) are related to the FFD y(t) throughthe linear system equations (167), we can describe the power spectra Sȳȳ (ω)
85and S∆ȳ∆ȳ (ω) through the power spectrum Syy (ω) of the FFD y(t), which
we introduced in equation (147). For this purpose, let us introduce system
functions A(ω) and B(ω) of the impulse responses a(t) and b(t), respectively:A(ω) =

Z∞−∞

B(ω) =

Z∞T

1a(t) e−iωt dt =Tb(t) e−iωt dt =−∞Z2e−iωt dt =

− T2


 Z0
1 

)
sin( ωT
2
,
ωT
( 2 )

e−iωt dt −

T 


−T

ZT

e−iωt dt





= 2i




0

)
sin2 ( ωT
2
.
ωT
( 2 )
(172)

From general properties of linear systems, we haveSȳȳ (ω) = Syy (ω) | A(ω) |2 ,andS∆ȳ∆ȳ (ω) = Syy (ω) | B(ω) |2 .
Therefore, equations (170) and (171) are reduced to
Squared system function of Allan variance
1.2

0.8

1
sin4(πντ)2(πντ)2sin2(πντ)(πντ)2

Squared system function of true variance
1

0.6
0.4
0.2
0

0.8
0.6
0.4
0.2

0

0.5

1

1.5

2

2.5

3

3.5

0

4

0

0.5

1

1.5

ντ

2
ντ

2.5

3

3.5

4

Figure 33: Squared system functions for the TVAR, | A(ω) |2 , (left) andthe AVAR, | B(ω) |2 , (right). Horizontal axes show the product ν T of thefrequency and the time interval.
1
I (T ) =2π2σy2 (T ) =

1
2π

Z∞

−∞
Z∞

Syy (ω))sin2 ( ωT2
ωT 2 dω,
( 2 )

(173)

sin4 ( ωT)2
ωT 2 dω.
( 2 )

(174)

Syy (ω) 2−∞

86

These are the equations which relate the TVAR and the AVAR to the power
spectrum Syy (ω) of the FFD y(t).
Figure 33 shows squared system functions | A(ω) |2 and | B(ω) |2 for the
TVAR and the AVAR, respectively. It is evident from this Figure that the
low frequency noise is effectively suppressed in the AVAR.
1.2.7

Time–Interval Dependence of Allan Variance

We can calculate dependence of the TVAR and the AVAR on the time interval
T according to equations (173) and (174) for each power–law index α of the
FFD spectrum Syy (ω) ∝ ω α given in Table 4. Results are shown in Table 5,
where ωh in the white phase (α = 2) and the flicker phase (α = 1) components
is a cut–off frequency at the high frequency side of the FFD spectrum.
Noise TypeSyy (ω)Sφφ (ω)White phaseH2 ω 2

ω02 H2Flicker phaseH1 ω 1

ω02 H1 ω −1

White frequency

H0

ω02 H0 ω −2

Flicker frequency

H−1 ω −1

ω02 H−1 ω −3

Random walk of frequency

H−2 ω −2

ω02 H−2 ω −4

σy2 (T )I 2 (T )3H2 ωhπT 2
3H1 ln(ωh T )
πT 2
H0
T
(2 ln 2) H−1
π
T H−2
3

2H2 ωh
πT 2

—
H0
T

—
—

Table 5: Time–interval dependence of the AVAR, σy2 (T ), and the TVAR,
I 2 (T ), for each power–law component of the FFD spectrum. ωh is a cut–off
frequency at the high frequency side of a power–law spectrum Syy (ω).
Figure 34 schematically shows the time–interval dependence of the “Allanstandard deviation (ASD)”, σy (T ), which is defined as the square root of the
AVAR. It is clear from this Figure that, having measured ASD as a functionof the time interval between samples T , we can readily distinguish noise
types, except for the case of the white phase and the flicker phase. This
demonstrates the real usefulness of the ASD (or the AVAR) as a measure of
the frequency stability.
Figure 35 shows characteristic performance of various frequency standards
in terms of the ASD versus the time interval. ASD’s of frequency standards
are approximated by different noise types in different time–interval ranges.
The active hydrogen maser, which shows ASD ≤ 10−15 in its flicker frequency
regime, exhibits the highest frequency stability in the time–interval range of
1,000 to 10,000 seconds.
87

Figure 34: A shematic view of the time–interval T dependence of the Allanstandard deviation σy (T ) in log–scales. Each noise type shows its own gradient in this diagram, though lines of the white phase and the flicker phase
are hardly distinguished by the difference of their gradients.

Figure 35: Performance of various frequency standards in terms of the ASD—
the time interval between samples. QZ: Quartz, RB: Rubidium, CS: Cesium,
HM: Active Hydrogen Maser. (Figure courtesy of HEWLETT PACKARD,
Application Note 1289).
881.2.8

True Variance and Allan Variance through Autocorrelationof Phase Noise

We can describe the TVAR and the AVAR also through the autocorrelation
Rφφ (T ) of the phase noise φ(t). In fact, for the TVAR, equations (154) and(158) yieldI 2 (T ) = hȳ 2 [k]i =2 [Rφφ (0) − Rφφ (T )]h[φ(tk + T ) − φ(tk )]2 i=.
2(ω0 T )
(ω0 T )2

(175)

Also, for the AVAR, equations (162) and (163) lead toh(∆ȳ[k])2 ih[φ(tk + 2T ) − 2 φ(tk + T ) + φ(tk )]2 i=
2
2 (ω0 T )2
3 Rφφ (0) − 4 Rφφ (T ) + Rφφ (2 T ).
=(ω0 T )2σy2 (T ) =(176)

Now, from equation (175), we have,1[R(0) − R(2 T )]2
I 2 (T ) − I 2 (2T ) =(ω0 T )2
3 Rφφ (0) − 4 Rφφ (T ) + Rφφ (2 T ).
=
2 (ω0 T )22 [R(0) − R(T )] −(177)

Therefore, in cases when the TVAR I 2 (T ) does not diverge, we have thefollowing relationship between the AVAR and the TVAR:
σy2 (T ) = 2 [I 2 (T ) − I 2 (2T )].
1.2.9

(178)

Coherence Function

As we saw earlier, the phase noise in the fringe phase causes the serious
problem, the coherence loss, for the sensitivity of VLBI.
The coherence loss due to the phase noise φ(t), after integrating the correlator output for duration T , can be estimated by introducing the “coherencefunction” (Rogers and Moran, 1981), which is defined byT

1 Z iφ(t)C(T ) =
e
dt .
T

(179)

0

Magnitude of the squared coherence function is represented by its dispersion:
T T

1 Z Z i[φ(t)−φ(t0 )]hC (T )i = 2hei dt dt0 .
T
2

0 0

89

(180)

Now, if we assume the Gaussian distribution of the phase noise differenceΦ = φ(t) − φ(t0 ):
Φ2
1
(181)
f (Φ) = √ e− 2σ2 ,
σ 2π
where
σ 2 = hΦ2 i = h[φ(t) − φ(t0 )]2 i,then, using the formulaZ∞2e−x −iax dx =√a2

π e− 4 ,−∞we obtainiΦ

he i =

Z∞

hΦ2 i

σ2

f (Φ) eiΦ dΦ = e− 2 = e− 2 .

−∞Hence, we reduce equation (180) to
T T

T T

1 Z Z − h[φ(t)−φ(t0 )]2 i
1 Z Z − D(t, t0 )2ee 2 dt dt0 ,hC (T )i = 2
dt dt0 = 2
T
T
2

0 0

(182)

0 0

where D(t, t0 ) is the so–called “temporal structure function” defined by
D(t, t0 ) = h[φ(t) − φ(t0 )]2 i.

(183)

Under the assumption of the stationary random phase noise, we haveD(t, t0 ) = D(τ ),with τ = t − t0 . Then we can reduce the double integral in equation(182)√
to a single integral. In fact, noting that dz = √12 dτ and L = 2 (T − τ ) in
Figure 36, we obtain
T

T

D(τ )
2 Z
2Z
− 2
hC (T )i = 2 (T − τ ) e
dτ =
T
T

2

0

0

τ
1−
T



D(τ )

e− 2 dτ.

(184)

According to equation (154), the mean FFD ȳ[k] with the sample interval
τ isφ(tk + τ ) − φ(tk ).
ȳ[k] =
ω0 τ
Therefore, we can describe the structure function D(τ ) through the TVAR
I 2 (τ ) asD(τ ) = ω02 τ 2 hȳ 2 i = ω02 τ 2 I 2 (τ ).
(185)
90

t’
T
L

T-τ

dz

τ
t

0
T

Figure 36: Geometry of integration in equation (182).
This enable us to describe the dispersion of the coherence function throughthe TVAR:
2hC (T )i =
T

ZT 

τ
1−
T



e−

2
=
T

ZT 

τ
1−
T



−

2

0

0

e

2 τ 2 I 2 (τ )
ω0
2

2 τ2
ω0
4π

R∞

−∞

dτ

Syy (ω)sin2 ( ωτ )2dω( ωτ )2
2

dτ.

(186)

Thus, when a functional form of the TVAR I 2 (τ ) is given, we can calculate the dispersion of the coherence function hC 2 (T )i . Then, we can estimate
the “coherence time τc ”, i.e. the interval during which we can more orless coherently integrate our signal. As the criterion for the coherence time
τc , we usually adopt the time interval which gives hC 2 (τc )i ∼ 0.85.
1.2.10Approximate Estimation of Coherence Time

Precise estimation of the coherence time using equation (186) is sometimes
impractical when the estimation of the TVAR with equation (159) diverses.
For practical purposes, the coherence time is usually estimated by a simpler way.
Phase noise accumulated during a time T is approximately given by '
ω0 σy T , where σy is “some standard deviation” of the FFD, for which weusually adopt the ASD.
It is obvious that we will not obtain any meaningful correlation result if
the phase noise varies more than 2π, during an integration time T . It is the91usual practice to require that the accumulated phase noise must not exceed
1 radian:
ω0 σy (T )T ≤ 1 radian.
(187)(Note that σy is, in general, a function of the time interval T .) Thus, thecoherence time τc is estimated as a time which satisfies
ω0 σy (τc )τc ' 1.

(188)

If we observe at ν = 8 GHz, and require that τc = 1000 sec, say, we needa frequency stability better than1
' 2 × 10−14 ,93
2π × 8 × 10 × 10
which is sometimes described as the “stability of a clock which would
deviate by 1 sec in 5 × 1013 sec, or 1.6 million years!”σy (1000 sec) 'The Active Hydrogen Maser Frequency Standard, which has the frequency stability of 10−16 < σy < 10−14 at time scales around 1000 sec (seeFigure 35), fulfilled the requirement, and is widely used in the world VLBI
observations.
1.2.11

Estimation of Time–Averaged Phase Noise

In some applications of VLBI obserrvation, for example in VLBI astrometry,it is meaningful to theoretically estimate the phase noise expected after time–
averaging of the correlator output for a duration of time T .
If we have a time series of measured fringe phase with high enough signal–
to–noise ratio, the time–averaged phase noise φ̄ is given by a simple model:
1
φ̄ =
T

ZT

φ(t) dt,(189)

0

where φ(t) is the phase noise at time t. This model is good enough for thecase of the “vector averaging” of the correlator output, such as shown inequation (179) for the coherence function, as long as the phase noise is kept
well smaller than 1 radian. In fact, if φ(t)  1, we have
1
T

ZT
0

e

iφ(t)

1
dt ∼
=
T

ZT

(1 + i φ(t)) dt ∼= eiφ̄ .

0Now, we will derive formulae for the dispersion of the residual phase noisefrom its time average:
1∆φ = φ(t) − φ̄ = φ(t) −
T
92

Z T
0

φ(t0 )dt0 ,(190)

and that of the time-averaged phase itself, given by equation (189).
The dispersions are
#2

"

1ZT
1ZT
φ(t) −
φ(t0 )dt0 dt
T 0
T 0
1 ZT
2 ZT
(T − τ )[Rφφ (0) − Rφφ (τ )]dτ = 2(T − τ ) ω02 τ 2 I 2 (τ )dτ.
= 2
T 0
T 0
(191)

2
σ∆φ
(T ) =

σφ̄2 (T ) =

"

1
T

Z T
0

φ(t)dt

#2 

=

2
T2

Z T
0



(T − τ )Rφφ (τ )dτ.

(192)

Therefore, we can calculate these dispersions if we know the TVAR I 2 (τ )and the autocorrelation of the phase noise Rφφ (τ ). Note that equation (191)
describes the accumulated phase noise around its time average. In particular,
equation (186) can be reduced to2hC (T )i =
T

ZT 

1−

τ
T



e−

2
∼
=
T

ZT 

τ
1−
T



[1 −

2

0

0

ω 2 τ 2 I 2 (τ )02dτω02 τ 2 I 2 (τ )2]dτ = 1 − σ∆φ
(T ),
2

(193)

as long as ω02 τ 2 I 2 (τ ) / 2  1.

In the simplest case of the white phase noise, we haveI 2 (τ ) =2 H 2 ωh,π τ2and σy2 (τ ) =3 H 2 ωh
,
π τ2

(194)

for the TVAR and the AVAR, andSφφ (ω) = ω02 H2 ,1
and Rφφ (τ ) =2πZ∞

Sφφ (ω) eiωτ dω = ω02 H2 δ(τ ),

−∞

(195)
for the power spectrum and the autocorrelation of the phase noise. Then,equations (191) and (192) yield2σ∆φ(T ) =

1 ZT
H2 ω02 ωh
2 2 2
(T
−
τ
)
ω
τ
I
(τ
)dτ
=
,
0
T2 0
π

(196)

for the dispersion of the residual phase noise, and
σφ̄2 (T ) =

2
T2

Z T
0

(T − τ )Rφφ (τ )dτ =93

H2 ω02,T(197)

for the dispersion of the time–averaged phase noise.
It is evident that, for the white phase noise, the dispersion of the residual2
phase noise σ∆φ(T ) does not depend on the averaging time T . This dispersion
is readily estimated, if we have measured AVAR σy2 (τ ), from which we caneasily extract the coefficient H2 ωh .
On the other hand, the dispersion of the time–averaged phase noise σφ̄2 (T )decreases with increasing averaging time as ∝ 1/T , just like the thermal
noise. In this case, only H2 figures as the unknown coefficient which is notdirectly available from the measured AVAR alone.
We can estimate the cut–off frequency ωh , and then H2 , in the following
way. Suppose that actual measurements of the phase noise are performed
with a time interval τm . If we regard the measured phase noise value asa kind of the running mean for the duration τm , the high frequency cut–off
should be ωh ∼= 1/τm . Then, we can estimate σφ̄2 (T ), using this ωh and H2 ωhderived from the measured AVAR.

1.3

Time Synchronization

In order to find the white fringe within the coherence interval, the clocksof the element antennas must be synchronized with accuracy τsync , whichshould be better than the coherence interval 2/B, where B is the recorded
bandwidth, as we saw in Chapter 3.
Therefore, even in the early days of VLBI observations, with the typicalbandwidth of B = 2 MHz, we needed high time–synchronization accuracy:
τsync < 1 µsec for B = 2MHz,
which was not easily available from time transfer systems at that time using
surface waves.
Nowadays, the requirement is much more severe since we use the observing bandwidth with B = 256 MHz or wider. Then, we need
τsync < 7.8 nsec(!) for B = 256MHz.
In actuality, multi–lag correlators, which we will discuss later, can significantly ease this requirement. But, anyway, we need highly accurate timesynchronization better than 100 nsec.
Fortunately, now GPS (Global Positioning System) Satellites are capable of providing time synchronizations at a few tens nsec level. Therefore,right now we do not have any essential problem in the time synchronization
technology.
Once VLBI fringe is successfully detected, VLBI itself serves as the best
time synchronization device, with 1 nsec level accuracy.
94

1.4

Recording System

Even in the very early stages, VLBI required a very high data rate andcapacity for the recording system. For example, if we wish to record digitized
data, with bandwidth B = 2 MHz, sampled at Nyquist rate (which equals
2B samples/sec) with 1 bit (or 2–level) quantization, for a time duration of
400 sec, we need a data recording rate of 4 Mbit/sec, and a data capacity ofat least 1.6 Gbit.
Nowadays, new recording sytems, such as the VERA system, allow 1 hour
recording at a 1 Gbit/sec rate, per volume of magnetic tape. Therefore, sucha tape records 3.6 Tbit of data.

Figure 37: VLBI systems (from Thompson, Moran, and Swenson, 2001).
Generations of “VLBI systems” were marked by the development of thedigital data recording technology, as shown in Figure 37.
The Communications Research Laboratory (CRL, now NICT: National
Institute of Information and Communications Technology) has developed theMark II–compatible K-1, the experimental real–time VLBI system K-2, theMark III–compatible K-3, and cassette–based K-4 systems, by its own efforts(Takahashi, et al., 2000).
Figures 38, 39, 40 and 41 show the VLBA recorder, S2 recorder, K-4
recorder, and the tape handler for K-4 recorder. Figure 42 shows the S2,
K-4, and VLBA (Mark IV) tapes.

95

Figure 38: VLBA recorder (MIT Haystack Observatory/NRAO, USA).

96

Figure 39: S2 recorder (Center for Research in Earth and Space Technology,
Canada).

97

Figure 40: K-4 recorder (NICT/CRL, Japan).

98

Figure 41: Automated tape handler for K-4 (NICT/CRL, Japan).

99

Figure 42: Magnetic tapes used in VLBI: S2 (top left), K-4 (left), and
VLBA/Mark IV (right).

100

2

Overview of the VLBI System

2.1

MK–3 As a Prototype of Modern VLBI Systems

Figure 43 shows a schematic view of the Mark III / K-3 VLBI system, often
called “MK–3”, which was initially developed by the MIT Haystack Observatory with sponsorship from NASA, and by the CRL (now NICT), in thelate 1970’s. Although this system has now been almost replaced by newer
systems, such as the Mark IV, Mark V, VLBA, K-4, K-5 and VSOP, many
important elements of these modern VLBI systems were implemented for the
first time in the MK–3 system, and have been further developed in the latest
systems. Therefore, we briefly describe here the basic components, and theirfunctions, of the MK–3 system.

P & D Antenna
Unit

S and X Low
Noise Amplifiers

Local Oscillator
S X
P&D
GroundUnit
H-Maser Frequency
Standard

Antenna
Control
Unit

Downconverter
Receiver RoomIF Frequency
Distributer VLBI Terminal
Video Converter
Formatter

ControlComputer

Data Recorder

Control Building

Figure 43: A schematic view of the Mark III / K-3 VLBI system.
The MK–3 system was originally developed mainly for realizing high–101

precision VLBI geodesy, although the system has been extensively used for
astrophysical radio source imaging observations as well. In order to accurately estimate the group delay, the system adopted a multi–frequency channel design, together with high–speed recording technology; these features
have been retained in most of the latest systems.
2.1.1

Dual–Frequency Reception

In geodetic applications, a VLBI antenna receives radio waves from distantextragalactic sources usually at two frequencies simultaneously, in order tocorrect the effects of frequency–dependent propagation delay in the ionosphere. The two frequencies, S–band (2 GHz) and X–band (8 GHz), aremost widely used in global geodetic VLBI observations. The RF frequency
bands in the MK–3 system typically covered 200 MHz for S–band (e.g., from2120 to 2320 MHz), and 420 MHz for X–band (e.g., from 8180 to 8600 MHz).
Multi–frequency coaxial feed horns are widely used for dual–frequency reception. Other systems, including FSS (Frequency Selective Surfaces) and spiralarrays (see, e.g., Figure 44), are also used.

Figure 44: S/X spiral–array feed system, developed by Hosei University,
Japan, is used for dual–frequency reception in VERA.

2.1.2

First Frequency Conversion

In the receiver room, usually built in the antenna structure, the RF signals
are amplified by the S–band and X–band low–noise amplifiers, and down102

converted to, for example, 100 – 300 MHz (S–band) and 100 – 520 MHz(X–band), respectively, with the local oscillator signals at frequencies of 2020
MHz (S–band) and 8080 MHz (X–band), which are generated from a reference signal provided by the Hydrogen Maser Frequency Standard.
2.1.3

Transmission to Control Building

Then, the IF signals are fed to a so–called MK–3 VLBI terminal rack in
a control building, via transmission cables which are usually laid in underground ducts, in order to minimize phase fluctuations due to cable lengthvariations, induced by temperature changes.
2.1.4

Intermediate Frequency Distributer
IF signal from receiver.

νIF Distributer
Makes n identical copies.

Video Converter
Selects specified frequency ranges
and converts them to video bands.

Formatter
Digitizes, formats and adds time codes.

to VLBI RecorderFigure 45: Data flow in the MK–3 VLBI terminal.
In the VLBI terminal rack, the S–band and X–band IF signals are copied103

over to 6 (S–band) and 8 (X–band) identical signals, respectively, by theso–called “Intermediate Frequency Distributer” unit (Figure 45).
2.1.5

Baseband Conversion

The 6 and 8 identical IF signals are then sent to 14 special units consisting of
frequency downconverters, called “Video Converters” (VC’s) (or “Baseband
Converters”, BBC’s) located in the same rack. Each Video Converter selects
a specified frequency range in the IF band, and downconverts this range intoa video band, typically from 0 to 2 MHz. For this purpose, a frequency synthesizer, which is built in the Video Converter, generates an LO frequency
at the edge of the selected IF range, using a reference signal provided by the
Hydrogen Maser Frequency Standard. The Video Converter is equipped witha sideband–rejection mixer, which is capable of converting both lower andupper sidebands around the selected LO frequency, separately. Thus, in total
up to 28 (14 × 2 SB’s) baseband frequency channels are generated by the
14 Video Converters. Usually, of course, the baseband channels are selectedin such a way that they correspond to different RF frequency ranges. Table
6 and Figure 46 show an example of the LO frequency distribution used ina recent geodetic VLBI observation. This apparently strange distribution of
the data among the baseband channels, which are to be subsequently digitized and recorded, are spread over fairly wide ranges of the RF bands; Thisis specifically designed for better estimation of the group delay observables
via a technique called “bandwidth synthesis”, which will be discussed later.

8200

8250

2220

2240

8300

2260

8350
8400
8450
RF frequency (MHz)

2280
2300
RF frequency (MHz)

8500

2320

8550

8600

2340

2360

Figure 46: Distributions of baseband channels within the RF frequency
ranges, in a sample geodetic VLBI observation, as used in Table 6. Top:X–band, and bottom: S–band.

104

Ch#

band

LO freequency (MHz)bandwidth (MHz)1

X

8209.99

4.00

2

X

8219.99

4.00

3

X

8249.99

4.00

4

X

8309.99

4.00

5

X

8419.99

4.00

6

X

8499.99

4.00

7

X

8549.99

4.00

8

X

8569.99

4.00

9

S

2239.99

4.00

10

S

2244.99

4.00

11

S

2259.99

4.00

12

S

2289.99

4.00

13

S

2319.99

4.00

14

S

2329.99

4.00

15

S

2344.99

4.00

16

S

2354.99

4.00

Table 6: An example of the distribution of LO frequencies in baseband channels. Note that this example is taken from a recent geodetic VLBI observation, which uses the K-4 VLBI system. The K-4 is newer than the MK–3,
and has 16 Video Converters (therefore, 16 baseband channels), each with 4
MHz bandwidth.

105

2.1.6

Formatter

The baseband signals from the 14 Video Converters (possibly containing
both USB and LSB channels) are then fed to 14 “Formatter” units in thesame rack, each of which converts the analog baseband signal into a digital
signal, using a high–speed sampler, with a one–bit (or two–level) quantization
scheme, and Nyquist sampling rate (i.e. 2BV C samples/sec, where BV C ≤ 2
MHz is the bandwidth of each video band channel). At the same time, theFormatter periodically generates time mark codes, using the reference signal
and the clock pulse provided by the Hydrogen Maser Frequency Standard,
and inserts them into the digitized data. These time mark codes play basic
and important roles in the later digital correlation processing.
2.1.7

Data Recorder

Up to 28 channels of the digitized, formatted data are sent to a high–speed
data recording device, the Honeywell M–96 open–reel digital recorder, in thecase of the MK–3 system and its successors (see Figure 38). The formatted
multi–channel data are recorded in parallel tracks.
From the above discussion, we see, that the maximum recording rate in
the MK–3 system is given by:
1 bit/sample × 2 × 2 M sample/sec (Nyquist rate) ×
2 sidebands × 14 Video Converters = 112 Mbit/sec.
2.1.8

Phase and Delay Calibration

The 14 Video Converters (VC’s) add their own arbitrary initial phases to theconverted video band signals. It is necessary, for the accurate estimation of
the group delay, to calibrate the phase offsets among the video bands due tothe VC initial phases. For that purpose, a special device called the “Phaseand Delay Calibrator” is used (Figure 47). This calibrator system, which consists of a Ground Unit in the VLBI terminal rack, and an Antenna Unit in
the receiver room, generates so-called “comb–tone” signal, which is nothing
but the comb function given in equations (19) and (23), i.e. equally spaced intime delta functions (“pulse series”). As we saw before, Fourier transform of
the comb function in the time domain is the comb function in the frequency
domain. Therefore, in the frequency domain, the comb–tone signal consistsof a large number of sine waves at equally spaced frequencies. Moreover, thesine waves have regularly aligned phases, based on the reference signal fromthe Hydrogen Maser Frequency Standard. The Antenna unit generates thecomb–tones and adds them to the observed signal in the RF band. In the106actual observations, the LO frequencies of the Video Converters are selectedin such a way that at least one comb–tone signal falls into each of the video
bands (e.g., at a frequency of 10 kHz in the MK–3 system). The comb–tone
signals pass through the same units as the observed data do, and, therefore,
are affected by the same phase offsets due to the VC initial phases, as theobserved data are. Then, the data containing the comb–tones are digitized
and recorded. These comb–tones are later detected in the correlation processing, and the phases of the comb–tones in the different video bands arecompared with each other, in order to estimate, and then remove from theobserved data, the unknown phase offsets due to the VC initial phases.
Antenna feedΣn Bn cos [nω(t-τ)]

A cos [ω(t-τ)]

Comb pulses in time

Antenna Unitt

Comb tones go throughthe same path as theradio source signal.

Round-trip
measurementof cable delay τ

IF Frequency
Distributer

Fourier transformωComb tones in spectrum

Video Converter
A cos (ωt)
Ground Unit

Formatter

H-MaserFrequency
Standard

Data Recorder

Comb tones at 10 kHz of
every video band channel

10 kHztCorrelation ProcessorωPhase of the 10 kHz tone is detected bycorrelating with a standard 10 kHz signal.

Figure 47: Phase and delay calibration system in the MK–3.

107

2.1.9Hydrogen Maser Frequency Standard

The Hydrogen Maser Frequency Standard provides very high stability reference signals, usually at 5 MHz or 10 MHz, to the first downconverters in thereceiver room, the Video Converters, the Formatters, and to the Phase andDelay Calibrator. It also provides clock pulses, usually at 1 PPS (pulse persecond), to many devices, to guarantee their synchronous operation. The Hydrogen Maser Frequency Standard is usually placed in a special magnetically–shielded and temperature–controlled room, to shied against external disturbances as much as possible. Figure 48 shows a Hydrogen Maser Frequency
Standard installed in the VERA Mizusawa station.

Figure 48: A Hydrogen Maser Frequency Standard in the VERA Mizusawa
station.

2.1.10Automated Operation by Control Computer

All devices in a VLBI observing station are usually designed to be fully controllable by a single control computer (Figure 49). All commands for setting
LO frequencies, sampler modes, recording rates, etc., are remotely controlled.
Also, commands for antenna operations, measurements of system noise temperature and meteorological parameters, and so on, are issued from the control computer. In a VLBI observation, a so–called VLBI schedule file is given
to the control computer, where all setting information, coordinates of the observed radio sources, and station locations are listed. Also, the detailed time
sequence of the observing events, such as repointing the antenna to the next108

source, tape start, tape stop, and so on, are given in Universal Times (UT).

P & D Antenna
Unit

S and X Low
Noise Amplifiers

Local Oscillator
S X
P&D
GroundUnit
H-Maser Frequency
Standard

Antenna
Control
Unit

Downconverter
Receiver RoomIF Frequency
Distributer VLBI Terminal
Video Converter
Formatter
Data Recorder

ControlComputer
Log
FileMark IIIDB

Schedule
File
Observation sequence
Source information
Antenna information

Control Building

Field System
Figure 49: A typical automated VLBI operation system.
The control computer automatically conducts all steps of the observation
by issuing commands to the antenna control unit, VLBI terminal, and data
recorder, according to the schedule file. As a result, the only remaining task,
to be done by an operator in a normal VLBI observation, is to change therecording tapes once an hour or, at most a few times a day, depending on thesystem design. The most widely used software for automated control of the
VLBI equipment, which was developed in the NASA Goddard Space Flight
Center (NASA/GSFC), is called the “VLBI Field System”.

2.2

Modern VLBI Systems

2.2.1

New Recording and Fiber–Link Systems

The major successors of the MK–3 system are listed in Figure 37. They are
the VLBA and Mark IV sytems based on the advanced head controls of the109

Honeywell M–96 recorder (Figure 38), the S2 system using 8 video cassetterecorders in parallel (Figure 39), and the K-4 system based on an ID1 digital
video cassette recorder (Figure 40).
Nowadays, modern VLBI recording systems with 1 Gbps recording speed
have been developed and some of them are already in use (Table 7). They
are the Mark V system (Figure 50) using a hard disk array instead of themagnetic tapes, the K–5 system also based on the hard disk array (Figure51), the GBR-1000 system based on an HDTV video cassette recorder (Figure52), the VERA recording system based on a new ID1 digital video cassetterecorder (Figure 53), and the fiber–linked e–VLBI.

name

type

bitrateMark V

hard disk array1024 Mbps

K–5

hard disk array512 & 2048 Mbps

GBR-1000

HDTV

1024 Mbps

VERA

ID1

1024 Mbps

e–VLBI

real–time fiber link

> 2048 Mbps

Table 7: New generation recording and fiber–link systems.

Figure 50: Mark V recording system based on a new hard disk array design.

110Figure 51: K–5 VLBI terminal.

111

Figure 52: GBR-1000 1 Gbps VLBI system based on an HDTV video cassetterecorder.

112Figure 53: VERA recording system based on a new ID1 digital video cassetterecorder.
2.2.2

Digital Baseband Converters

Another progress is beeing made in several countries for implementing digital
BBC systems based on the digital filtering technology. The idea is first to

Figure 54: Data Aquisition System of the KVN adopting the digital filter
system for the baseband conversion.
digitize a wideband IF data using a high–speed sampler (analog–to–digitalconverter), and second to cut the digital data into baseband channels bymeans of a high–speed digital filter. Then the baseband data will no longer
suffer from irregular bandpass characteristics of analog BBC’s and arbitrary113phases added by the LO’s of the BBC’s. In fact, the digital filter can in principle remove the phase calibration systems discussed above, as far as the LOphase problem only is concerned. Also, the digital filter allows flexible organization of baseband channels meeting different scientific requirements. For
example, a single channel ultra–wideband data could be suited to imaging
very weak continuum sources, while 16 channel data are necessary for geodetic observations. Such a digital BBC system was developed and successfully
utilized in 4 stations of the VERA array.
Figure 54 shows the Data Aqusition System (DAS) of the KVN whichuses a digital BBC system for the baseband conversion. The DAS system
includes 4 high–speed samplers located in the receiver room. The 4 samplers
are equipped for 4 frequency bands to be simultaneously received for mm–
wave VLBI observations. Each sampler is capable of digitizing IF data of
512 MHz bandwidth using the higher–order sampling technique discussed
earlier with 1 Gsps speed in the 2–bit quantization mode (therefore, theoutput bit–rate is 2 Gbps). The digital data are transmitted to the control
building via optical fiber cables. Then the digital filter forms 1, 2, 4, 8, and16 baseband channels with bandwidths of 256, 128, 64, 32, and 16 or 8 MHz,correspondingly, out of the 512 MHz / 2–bit / 1 Gsps input, depending on
scientific purposes.
2.2.3

e–VLBI

Several groups in the world are now developing VLBI systems based on high–speed data transmission techniques via optical fibers. The idea is to replace
the data tapes by ultra–wide–band transmission cables. This technology
is now called “e–VLBI”. Figures 55 and 56 show the first successful EVN(European VLBI Network) e–VLBI observation conducted by three observatories, in April 2004, and a beautiful image of the gravitational lens object
B0218+357 obtained through this coordinated effort.
Of course, this technology by no means transforms VLBI into a connected–
element interferometer, since the frequency standards in different stations
must remain independent, in view of the essential technical difficulties totransmit reference signals over thousands of kilometers without significant
delays or phase fluctuations. Nevertheless, this technology will bring VLBI
much closer to the connected–element interferometer, in the sense that theobserved data could be correlated and analyzed in real–time, or almost inreal time (“near–real-time”). For example, pioneering real–time e–VLBI experiments conducted in the KSP (Key Stone Project, 1995–2001) of the CRL(now NICT) regularly yielded final geodetic results a few minutes after eachobservation. Moreover, the optical fiber cables offer even higher data trans114

Figure 55: First EVN e–VLBI observation (2004 April) using radio telescopes
at Jodrell Bank, UK; Westerbork, the Netherlands; and Onsala, Sweden. The
lower panels show beats of the fringe amplitudes produced by two closely
spaced sources, and an image map obtained during the observation (from
URL: http://www.jive.nl).

Figure 56: A close–up view of the image of a radio–loud gravitational lens
JVAS B0218+357, obtained in the first EVN e–VLBI observation (from URL:
http://www.jive.nl).
115

mission rates, and therefore higher sensitivity, than the magnetic tapes, or
hard disk arrays. Real–time VLBI experiments with 2.5 Gbps transmission
rate have been successfully conducted since 1998 (Figure 57). Even muchhigher transmission rates are expected in new connected–element interferometer arrays (for example, 96 Gbps / antenna is planned for the Expanded
VLA (EVLA), and for the Atacama Large Millimeter and submillimeter Array (ALMA)). High–speed correlators are now being developed to meet thesehigh data rates.
At the same time, less expensive and more widely accessible VLBI datatransmission, via broad–band Internet using the IP protocol, is also being
intensively studied, and has been successfully tested. This IP–based e–VLBI,
or IP–VLBI, seems to be a particularly promising technology which will
make VLBI observations much more user–friendly for many astronomers andgeophysicists around the world. The K–5 system, developed at the NICT, isdesigned to realize the IP–VLBI concept (Figure 58).
2.2.4

VLBI Standard Interface (VSI)

Another remarkable example of progress in modern VLBI is the definition of
the international “VLBI Standard Interface (VSI)” specifications(http://web.haystack.edu/vsi).
The world VLBI community has long suffered from incompatibily of thevarious different VLBI systems, shown above, which had been developed indifferent institutions. In order to cross–correlate data recorded on different
tapes by different VLBI systems A and B, say, one had to convert the formatof B’s data into A’s format, copy the converted data to A’s tape, and thencross–correlate them with A’s correlator, or vice versa.
The purpose of the VSI specifications is to make all VLBI systems in theworld compatible, provided only that they obey these standard specifications.
So far, a hardware specification called “VSI–H”, and a software specification called “VSI–S”, were worked out by an international working group of
VLBI specialists. Another specification for e–VLBI, to be called “VSI–E”,
is now under intensive discussion and development.
The VSI–H hardware specification defines a standard interface between
the “data–acquisition system (DAS)” (i.e. VLBI terminals) and the “data–transmission sytem (DTS)” (i.e. tapes or hard disks), as well as between the“data–transmission system” and the “data–processing system (DPS)” (i.e.
correlators). Figure 59 shows how the interface could be realized. Every
DAS, DTS, and DPS must have common connectors with well specified pin–assignments and data rates. The actual data transmission media (the handler
units for tapes or hard–disks) are supposed to be equipped with a “data116Figure 57: 2.5 Gbps fiber–linked real–time VLBI experiment “OLIVE–GALAXY”, Japan (1998, December).

117

Figure 58: IP–based e–VLBI developed for the K–5 system.

118

Figure 59: VSI-H functional diagram (http://web.haystack.edu/vsi).

119

input module (DIM)”, and a “data output module (DOM)”, which interface
the common connectors, and the connectors of the handler hardware. The
actual data formats in the individual data transmission media are arbitrary,provided that the DIM accepts the input data stream (specified by VSI-H),
and the DOM yields the output data stream (also specified by VSI-H).
The VSI–S software specification defines protocols for handling VSI–H–
compliant equipment.
A number of VSI–H–compliant VLBI systems have been developed invarious countries, and successfully cross–correlated with each other.

3

Difficulties in Ground–Based VLBIUNFINISHED.

4

Correlation Processing in VLBIUNFINISHED.

5

Observables of VLBIUNFINISHED.

120

References
Allan, D.W., 1966, Statistics of Atomic Frequency Standards,
Proc. IEEE, 54, 221–230.
Hagen, J.B., and Farley, D.T., 1973, Digital Correlation Techniques in Radio Science, Radio Sci., 8, 775–784.
Iguchi, S., and Kawaguchi, N., 2002, Higher–Order Sampling,
Over Sampling and Digital Filtering Techniques for Radio Interferometry, IEICE, Trans. Commun., E85-B, 1806–1816.
Meijering, E., 2002, A Chronology of Interpolation: From Ancient
Astronomy to Modern Signal and Image Processing, Proc,
IEEE, 90, 319–342.
Papoulis, A., 1984, Probability, Random Variables, and Stochastic Processes, 2nd Edition, McGraw-Hill (New York).
Price, R., 1958, A Useful Theorem for Nonlinear Devices Having
Gaussian Inputs, IRE Trans. Inf. Theory, IT-4, 69–72.
Rogers, A.E.E., and Moran, J.M., 1981, Coherence Limits for
Very–Long–Baseline Interferometry, IEEE Trans. Instrum.
Meas., IM–30, 283–286.
Shannon, C.E., 1949, Communication in the Presence of Noise,
Proc. IEEE, 37, 10–21.
Thompson, A.R., Moran, J.M., and Swenson, G.W., Jr., 2001,
Interferometry and Synthesis in Radio Astronomy, John Wiley& Sons (New York).
Takahashi, F., Kondo, T., Takahashi, Y., and Koyama, Y., 2000,
Very Long Baseline Interferometer, Ohmsha (Tokyo).
Van Vleck, J.H., and Middelton, D., 1966, The Spectrum of
Clipped Noise, Proc. IEEE, 54, 2–19.

121