# NAG C Library Function Document

## 1Purpose

nag_tsa_spectrum_bivar (g13cdc) calculates the smoothed sample cross spectrum of a bivariate time series using spectral smoothing by the trapezium frequency (Daniell) window.

## 2Specification

 #include #include
 void nag_tsa_spectrum_bivar (Integer nxy, NagMeanOrTrend mt_correction, double pxy, Integer mw, Integer is, double pw, Integer l, Integer kc, const double x[], const double y[], Complex **g, Integer *ng, NagError *fail)

## 3Description

The supplied time series may be mean and trend corrected and tapered as in the description of nag_tsa_spectrum_univar (g13cbc) before calculation of the unsmoothed sample cross-spectrum
 $f xy * ω = 1 2πn ∑ t=1 n y t expiωt × ∑ t=1 n x t exp -i ω t$
for frequency values ${\omega }_{j}=\frac{2\pi j}{K}$, $0\le {\omega }_{j}\le \pi$.
A correction is made for bias due to any tapering.
As in the description of nag_tsa_spectrum_univar (g13cbc) for univariate frequency window smoothing, the smoothed spectrum is returned at a subset of these frequencies,
 $ν l = 2πl L , l = 0 , 1 , … , L/2$
where [ ] denotes the integer part.
Its real part or co-spectrum $cf\left({\nu }_{l}\right)$, and imaginary part or quadrature spectrum $qf\left({\nu }_{l}\right)$ are defined by
 $f xy ν l = cf ν l + iqf ν l = ∑ ω k < π/M w ~ k f xy * ν l + ω k$
where the weights ${\stackrel{~}{w}}_{k}$ are similar to the weights ${w}_{k}$ defined for nag_tsa_spectrum_univar (g13cbc), but allow for an implicit alignment shift $S$ between the series:
 $w ~ k = w k exp -2 π iSk / L .$
It is recommended that $S$ is chosen as the lag $k$ at which the cross-covariances ${c}_{xy}\left(k\right)$ peak, so as to minimize bias.
If no smoothing is required, the integer $M$ which determines the frequency window width $\frac{2\pi }{M}$, should be set to $n$.
The bandwidth of the estimates will normally have been calculated in a previous call of nag_tsa_spectrum_univar (g13cbc) for estimating the univariate spectra of ${y}_{t}$ and ${x}_{t}$.
Bloomfield P (1976) Fourier Analysis of Time Series: An Introduction Wiley
Jenkins G M and Watts D G (1968) Spectral Analysis and its Applications Holden–Day

## 5Arguments

1:    $\mathbf{nxy}$IntegerInput
On entry: the length of the time series $x$ and $y$, $n$.
Constraint: ${\mathbf{nxy}}\ge 1$.
2:    $\mathbf{mt_correction}$NagMeanOrTrendInput
On entry: whether the data are to be initially mean or trend corrected. ${\mathbf{mt_correction}}=\mathrm{Nag_NoCorrection}$ for no correction, ${\mathbf{mt_correction}}=\mathrm{Nag_Mean}$ for mean correction, ${\mathbf{mt_correction}}=\mathrm{Nag_Trend}$ for trend correction.
Constraint: ${\mathbf{mt_correction}}=\mathrm{Nag_NoCorrection}$, $\mathrm{Nag_Mean}$ or $\mathrm{Nag_Trend}$.
3:    $\mathbf{pxy}$doubleInput
On entry: the proportion of the data (totalled over both ends) to be initially tapered by the split cosine bell taper.
A value of $0.0$ implies no tapering.
Constraint: $0.0\le {\mathbf{pxy}}\le 1.0$.
4:    $\mathbf{mw}$IntegerInput
On entry: the frequency width, $M$, of the smoothing window as $2\pi /M$.
A value of $n$ implies that no smoothing is to be carried out.
Constraint: $1\le {\mathbf{mw}}\le {\mathbf{nxy}}$.
5:    $\mathbf{is}$IntegerInput
On entry: the alignment shift, $S$, between the $x$ and $y$ series. If $x$ leads $y$, the shift is positive.
Constraint: $-{\mathbf{l}}<{\mathbf{is}}<{\mathbf{l}}$.
6:    $\mathbf{pw}$doubleInput
On entry: the shape argument, $p$, of the trapezium frequency window.
A value of $0.0$ gives a triangular window, and a value of $1.0$ a rectangular window.
If ${\mathbf{mw}}={\mathbf{nxy}}$ (i.e., no smoothing is carried out) then pw is not used.
Constraint: $0.0\le {\mathbf{pw}}\le 1.0$ if ${\mathbf{mw}}\ne {\mathbf{nxy}}$.
7:    $\mathbf{l}$IntegerInput
On entry: the frequency division, $L$, of smoothed cross spectral estimates as $2\pi /L$.
Constraint: ${\mathbf{l}}\ge 1$.
l must be a factor of kc (see below).
8:    $\mathbf{kc}$IntegerInput
On entry: the order of the fast Fourier transform (FFT) used to calculate the spectral estimates. kc should be a product of small primes such as ${2}^{m}$ where $m$ is the smallest integer such that ${2}^{m}\ge 2n$, provided $m\le 20$.
Constraints:
• ${\mathbf{kc}}\ge 2×{\mathbf{nxy}}$;
• kc must be a multiple of l. The largest prime factor of kc must not exceed $19$, and the total number of prime factors of kc, counting repetitions, must not exceed $20$. These two restrictions are imposed by the internal FFT algorithm used.
9:    $\mathbf{x}\left[{\mathbf{kc}}\right]$const doubleInput
On entry: the nxy data points of the $x$ series.
10:  $\mathbf{y}\left[{\mathbf{kc}}\right]$const doubleInput
On entry: the nxy data points of the $y$ series.
11:  $\mathbf{g}$Complex **Output
On exit: the complex vector which contains the ng cross spectral estimates in elements ${\mathbf{g}}\left[0\right]$ to ${\mathbf{g}}\left[{\mathbf{ng}}-1\right]$. The $y$ series leads the $x$ series.
The memory for this vector is allocated internally. If no memory is allocated to g (e.g., when an input error is detected) then g will be NULL on return. If repeated calls to this function are required then NAG_FREE should be used to free the memory in between calls.
12:  $\mathbf{ng}$Integer *Output
On exit: the number of spectral estimates, $\left[L/2\right]+1$, whose separate parts are held in g.
13:  $\mathbf{fail}$NagError *Input/Output
The NAG error argument (see Section 3.7 in How to Use the NAG Library and its Documentation).

## 6Error Indicators and Warnings

NE_2_INT_ARG_CONS
On entry, ${\mathbf{kc}}=〈\mathit{\text{value}}〉$ while ${\mathbf{l}}=〈\mathit{\text{value}}〉$. These arguments must satisfy kc%${\mathbf{l}}\ne 0$ when ${\mathbf{l}}>0$.
On entry, ${\mathbf{kc}}=〈\mathit{\text{value}}〉$ while ${\mathbf{nxy}}=〈\mathit{\text{value}}〉$. These arguments must satisfy ${\mathbf{kc}}\ge 2$*nxy when ${\mathbf{nxy}}>0$.
On entry, ${\mathbf{l}}=〈\mathit{\text{value}}〉$ while ${\mathbf{is}}=〈\mathit{\text{value}}〉$. These arguments must satisfy $‖{\mathbf{is}}‖<{\mathbf{l}}$ when ${\mathbf{l}}>0$.
NE_2_INT_ARG_GT
On entry, ${\mathbf{mw}}=〈\mathit{\text{value}}〉$ while ${\mathbf{nxy}}=〈\mathit{\text{value}}〉$. These arguments must satisfy ${\mathbf{mw}}\le {\mathbf{nxy}}$.
NE_ALLOC_FAIL
Dynamic memory allocation failed.
On entry, argument mt_correction had an illegal value.
NE_FACTOR_GT
At least one of the prime factors of kc is greater than $19$.
NE_INT_ARG_LT
On entry, ${\mathbf{l}}=〈\mathit{\text{value}}〉$.
Constraint: ${\mathbf{l}}\ge 1$.
On entry, ${\mathbf{mw}}=〈\mathit{\text{value}}〉$.
Constraint: ${\mathbf{mw}}\ge 1$.
On entry, ${\mathbf{nxy}}=〈\mathit{\text{value}}〉$.
Constraint: ${\mathbf{nxy}}\ge 1$.
NE_INTERNAL_ERROR
An internal error has occurred in this function. Check the function call and any array sizes. If the call is correct then please contact NAG for assistance.
NE_REAL_ARG_GT
On entry, pw must not be greater than 1.0: ${\mathbf{pw}}=〈\mathit{\text{value}}〉$.
On entry, pxy must not be greater than 1.0: ${\mathbf{pxy}}=〈\mathit{\text{value}}〉$.
NE_REAL_ARG_LT
On entry, pw must not be less than 0.0: ${\mathbf{pw}}=〈\mathit{\text{value}}〉$.
On entry, pxy must not be less than 0.0: ${\mathbf{pxy}}=〈\mathit{\text{value}}〉$.
NE_TOO_MANY_FACTORS
kc has more than 20 prime factors.

## 7Accuracy

The FFT is a numerically stable process, and any errors introduced during the computation will normally be insignificant compared with uncertainty in the data.

## 8Parallelism and Performance

nag_tsa_spectrum_bivar (g13cdc) is not threaded in any implementation.

nag_tsa_spectrum_bivar (g13cdc) carries out an FFT of length kc to calculate the sample cross spectrum. The time taken by the function for this is approximately proportional to ${\mathbf{kc}}×\mathrm{log}\left({\mathbf{kc}}\right)$ (but see function document nag_sum_fft_realherm_1d (c06pac) for further details).

## 10Example

The example program reads 2 time series of length $296$. It selects mean correction and a 10% tapering proportion. It selects a $2\pi /16$ frequency width of smoothing window, a window shape argument of $0.5$ and an alignment shift of $3$. It then calls nag_tsa_spectrum_bivar (g13cdc) to calculate the smoothed sample cross spectrum and prints the results.

### 10.1Program Text

Program Text (g13cdce.c)

### 10.2Program Data

Program Data (g13cdce.d)

### 10.3Program Results

Program Results (g13cdce.r)