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

# Syntax

C#
```public static void g13cd(
int nxy,
int mtxy,
double pxy,
int mw,
int ish,
double pw,
int l,
int kc,
double[] xg,
double[] yg,
out int ng,
out int ifail
)```
Visual Basic
```Public Shared Sub g13cd ( _
nxy As Integer, _
mtxy As Integer, _
pxy As Double, _
mw As Integer, _
ish As Integer, _
pw As Double, _
l As Integer, _
kc As Integer, _
xg As Double(), _
yg As Double(), _
<OutAttribute> ByRef ng As Integer, _
<OutAttribute> ByRef ifail As Integer _
)```
Visual C++
```public:
static void g13cd(
int nxy,
int mtxy,
double pxy,
int mw,
int ish,
double pw,
int l,
int kc,
array<double>^ xg,
array<double>^ yg,
[OutAttribute] int% ng,
[OutAttribute] int% ifail
)```
F#
```static member g13cd :
nxy : int *
mtxy : int *
pxy : float *
mw : int *
ish : int *
pw : float *
l : int *
kc : int *
xg : float[] *
yg : float[] *
ng : int byref *
ifail : int byref -> unit
```

#### Parameters

nxy
Type: System..::..Int32
On entry: $n$, the length of the time series $x$ and $y$.
Constraint: ${\mathbf{nxy}}\ge 1$.
mtxy
Type: System..::..Int32
On entry: whether the data is to be initially mean or trend corrected.
${\mathbf{mtxy}}=0$
For no correction.
${\mathbf{mtxy}}=1$
For mean correction.
${\mathbf{mtxy}}=2$
For trend correction.
Constraint: $0\le {\mathbf{mtxy}}\le 2$.
pxy
Type: System..::..Double
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$.
mw
Type: System..::..Int32
On entry: $M$, the frequency width of the smoothing window as $\frac{2\pi }{M}$.
A value of $n$ implies that no smoothing is to be carried out.
Constraint: $1\le {\mathbf{mw}}\le {\mathbf{nxy}}$.
ish
Type: System..::..Int32
On entry: $S$, the alignment shift between the $x$ and $y$ series. If $x$ leads $y$, the shift is positive.
Constraint: $-{\mathbf{l}}<{\mathbf{ish}}<{\mathbf{l}}$.
pw
Type: System..::..Double
On entry: $p$, the shape parameter 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: if ${\mathbf{mw}}\ne {\mathbf{nxy}}$, $0.0\le {\mathbf{pw}}\le 1.0$.
l
Type: System..::..Int32
On entry: $L$, the frequency division of smoothed cross spectral estimates as $\frac{2\pi }{L}$.
Constraints:
• ${\mathbf{l}}\ge 1$;
• l must be a factor of kc.
kc
Type: System..::..Int32
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.
xg
Type: array<System..::..Double>[]()[][]
An array of size [kc]
On entry: the nxy data points of the $x$ series.
On exit: the real parts of the ng cross spectral estimates in elements ${\mathbf{xg}}\left[0\right]$ to ${\mathbf{xg}}\left[{\mathbf{ng}}-1\right]$, and ${\mathbf{xg}}\left[{\mathbf{ng}}\right]$ to ${\mathbf{xg}}\left[{\mathbf{kc}}-1\right]$ contain $0.0$. The $y$ series leads the $x$ series.
yg
Type: array<System..::..Double>[]()[][]
An array of size [kc]
On entry: the nxy data points of the $y$ series.
On exit: the imaginary parts of the ng cross spectral estimates in elements ${\mathbf{yg}}\left[0\right]$ to ${\mathbf{yg}}\left[{\mathbf{ng}}-1\right]$, and ${\mathbf{yg}}\left[{\mathbf{ng}}\right]$ to ${\mathbf{yg}}\left[{\mathbf{kc}}-1\right]$ contain $0.0$. The $y$ series leads the $x$ series.
ng
Type: System..::..Int32%
On exit: the number of spectral estimates, $\left[L/2\right]+1$, whose separate parts are held in xg and yg.
ifail
Type: System..::..Int32%
On exit: ${\mathbf{ifail}}={0}$ unless the method detects an error or a warning has been flagged (see [Error Indicators and Warnings]).

# Description

The supplied time series may be mean and trend corrected and tapered as in the description of g13cb before calculation of the unsmoothed sample cross-spectrum
 $fxy*ω=12πn∑t=1nytexpiωt×∑t=1nxtexp-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 g13cb for univariate frequency window smoothing, the smoothed spectrum is returned as a subset of these frequencies,
 $νl=2πlL, 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
 $fxyνl=cfνl+iqfνl=∑ωk<πMw~kfxy*νl+ωk$
where the weights ${\stackrel{~}{w}}_{k}$ are similar to the weights ${w}_{k}$ defined for g13cb, but allow for an implicit alignment shift $S$ between the series:
 $w~k=wkexp-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 g13cb for estimating the univariate spectra of ${y}_{t}$ and ${x}_{t}$.

# References

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

# Error Indicators and Warnings

Errors or warnings detected by the method:
${\mathbf{ifail}}=1$
 On entry, ${\mathbf{nxy}}<1$, or ${\mathbf{mtxy}}<0$, or ${\mathbf{mtxy}}>2$, or ${\mathbf{pxy}}<0.0$, or ${\mathbf{pxy}}>1.0$, or ${\mathbf{mw}}<1$, or ${\mathbf{mw}}>{\mathbf{nxy}}$, or ${\mathbf{pw}}<0.0$ and ${\mathbf{mw}}\ne {\mathbf{nxy}}$, or ${\mathbf{pw}}>1.0$ and ${\mathbf{mw}}\ne {\mathbf{nxy}}$, or ${\mathbf{l}}<1$, or $\left|{\mathbf{ish}}\right|\ge {\mathbf{l}}$.
${\mathbf{ifail}}=2$
 On entry, ${\mathbf{kc}}<2×{\mathbf{nxy}}$, or kc is not a multiple of l, or kc has a prime factor exceeding $19$, or kc has more than $20$ prime factors, counting repetitions.
${\mathbf{ifail}}=3$
This indicates that a serious error has occurred. Check all array subscripts in calls to g13cd. Seek expert help.
${\mathbf{ifail}}=-9000$
An error occured, see message report.
${\mathbf{ifail}}=-8000$
Negative dimension for array $〈\mathit{\text{value}}〉$
${\mathbf{ifail}}=-6000$
Invalid Parameters $〈\mathit{\text{value}}〉$

# Accuracy

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

# Parallelism and Performance

None.

g13cd carries out an FFT of length kc to calculate the sample cross spectrum. The time taken by the method for this is approximately proportional to ${\mathbf{kc}}×\mathrm{log}\left({\mathbf{kc}}\right)$ (but see method document (C06PAF not in this release) for further details).

# Example

This example reads two 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 parameter of $0.5$ and an alignment shift of $3$. It then calls g13cd to calculate the smoothed sample cross spectrum and prints the results.

Example program (C#): g13cde.cs

Example program data: g13cde.d

Example program results: g13cde.r