Covariance Module

class kite.Covariance(scene, config=<kite.covariance.CovarianceConfig object>)[source]

Construct the variance-covariance matrix of quadtree subsampled data.

Variance and covariance estimates are used to construct the weighting matrix to be used later in an optimization.

Two different methods exist to propagate full-resolution data variances and covariances of kite.Scene.displacement to the covariance matrix of the subsampled dataset:

  1. The distance between kite.quadtree.QuadNode leaf focal points, kite.covariance.Covariance.matrix_focal defines the approximate covariance of the quadtree leaf pair.

  2. The _accurate_ propagation of covariances by taking the mean of every node pair pixel covariances. This process is computational very expensive and can take a few minutes. kite.covariance.Covariance.matrix_focal

Parameters
covariance_matrix
Getter

(Cached) Covariance matrix calculated from mean of all pixel pairs inside the node pairs (full and accurate propagation).

Type

numpy.ndarray, size (nleaves x nleaves)

covariance_matrix_focal
Getter

(Cached) Approximate Covariance matrix from quadtree leaf pair distance only. Fast, use for intermediate steps only and finallly use approach covariance_matrix.

Type

numpy.ndarray, size (nleaves x nleaves)

covariance_model
Covariance model parameters for

modelCovariance() retrieved from getCovarianceFunction.

Note

using this function implies several model fits: (1) fit of the spectrum and (2) the cosine transform. Not sure about the consequences, if this is useful and/or meaningful.

Getter

Get the parameters.

Type

tuple, a and b

covariance_model_rms
Getter

RMS missfit between covariance_model and getCovarianceFunction

Type

float

covariance_spatial
Getter

(Cached) Undocumented

covariance_spectral
Getter

(Cached) Covariance function estimated directly from the power spectrum of displacement noise patch using the cosine transform.

Type

tuple, numpy.ndarray (covariance, distance)

export_weight_matrix(filename)[source]
Export the full weight_matrix to an ASCII

file. The data can be loaded through numpy.loadtxt().

Parameters

filename (str) – path to export to

getCovariance()[source]

Calculate the covariance function

Returns

The covariance and distance

Return type

tuple

getLeafCovariance(leaf1, leaf2)[source]
Get the covariance between leaf1 and leaf2 from

distances.

Parameters
  • leaf1 (str of leaf.id or QuadNode) – Leaf one

  • leaf2 (str of leaf.id or QuadNode) – Leaf two

Returns

Covariance between leaf1 and leaf2

Return type

float

getLeafWeight(leaf, model='focal')[source]
Get the total weight of leaf, which is the summation of

all single pair weights of kite.Covariance.weight_matrix.

\[w_{x} = \sum_i W_{x,i}\]
Parameters
  • model (str) – Focal or full, default focal

  • leaf (QuadNode) – A leaf from Quadtree

Returns

Weight of the leaf

Return type

float

getQuadtreeNoise(rstate=None, gather=<function nanmedian>)[source]

Create noise for a Quadtree

Use getSyntheticNoise() to create data-driven noise on each quadtree leaf, summarized by

Parameters

gather – Function gathering leaf’s noise realisation, defaults to num.median.

Returns

Array of noise level at each quadtree leaf.

Return type

numpy.ndarray

getStructure(method=None)[source]

Get the structure function

Parameters

method (str (optional)) – Either spatial or spectral, if None the method is taken from config

Returns

(variance, distance)

Return type

tuple

noise_coord

Coordinates of the noise patch in local coordinates.

Setter

Set the noise coordinates

Getter

Get the noise coordinates

Type

numpy.ndarray, [llE, llN, sizeE, sizeN]

noise_patch_size_km2
Getter

Noise patch size in \(km^2\).

Type

float

nthreads
Number of threads (CPU cores) to use for full covariance

calculation

Setting nthreads to 0 uses all available cores (default).

Setter

Sets the number of threads

Type

int

plot
Getter

(Cached) Simple overview plot to summarize the covariance estimations.

plot_syntheticNoise
Getter

(Cached) Simple overview plot to summarize the covariance estimations.

selectNoiseNode()[source]

Choose noise node from quadtree the biggest QuadNode from Quadtree.

Returns

A quadnode with the least signal.

Return type

QuadNode

setConfig(config=None)[source]

Sets and updated the config of the instance

Parameters

config (CovarianceConfig, optional) – New config instance, defaults to configuration provided by parent Scene

setSamplingMethod(method)[source]

Set the sampling method

setSpatialBins(nbins)[source]

Set number of spatial bins

setSpatialPairs(npairs)[source]

Set number of random spatial pairs

structure_spatial
Getter

(Cached) Undocumented

structure_spectral
Getter

(Cached) Structure function derived from noise_patch :type: tuple, numpy.ndarray (structure_spectral, distance)

Adapted from http://clouds.eos.ubc.ca/~phil/courses/atsc500/docs/strfun.pdf

syntheticNoise(shape=(1024, 1024), dEdN=None, anisotropic=False, rstate=None)[source]

Create random synthetic noise from data noise power spectrum.

This function uses the power spectrum of the data noise (noise_data) (powerspecNoise()) to create synthetic noise, e.g. to use it for data pertubation in optinmizations. The default sampling distances are taken from kite.scene.Frame.dE and kite.scene.Frame.dN. They can be overwritten.

Parameters
  • shape (tuple, optional) – shape of the desired noise patch. Pixels in northing and easting (nE, nN), defaults to (1024, 1024).

  • dEdN – The sampling distance in easting, defaults to (kite.scene.Frame.dE, kite.scene.Frame.dN).

Returns

synthetic noise patch

Return type

numpy.ndarray

weight_matrix
Getter

(Cached) Weight matrix from full covariance \(cov^{-1}\).

Type

numpy.ndarray, size (nleaves x nleaves)

weight_matrix_L2
Getter

(Cached) Weight matrix from full covariance \(\sqrt{cov^{-1}}\).

Type

numpy.ndarray, size (nleaves x nleaves)

weight_matrix_focal
Getter

(Cached) Approximated weight matrix from fast focal method \(cov_{focal}^{-1}\).

Type

numpy.ndarray, size (nleaves x nleaves)

weight_vector
Getter

(Cached) Weight vector from full covariance \(cov^{-1}\).

Type

numpy.ndarray, size (nleaves)

weight_vector_focal
Getter

(Cached) Weight vector from fast focal method \(\sqrt{cov_{focal}^{-1}}\).

Type

numpy.ndarray, size (nleaves)

CovarianceConfig

class kite.covariance.CovarianceConfig(*args, **kwargs)[source]

Undocumented.

noise_coord

numpy.ndarray (pyrocko.guts_array.Array), optional

Noise patch coordinates and size,

model_coefficients

tuple of pyrocko.guts.Any objects, optional

Covariance model coefficients. Either two (exponential) or three (exponential and cosine term) coefficients.See also modelCovariance().

model_function

builtins.str (pyrocko.guts.StringChoice), default: 'exponential'

Covariance approximation function.

sampling_method

builtins.str (pyrocko.guts.StringChoice), default: 'spatial'

Method for estimating the covariance and structure function.

spatial_bins

int, default: 75

Number of distance bins for spatial covariance sampling.

spatial_pairs

int, default: 200000

Number of random pairs for spatial covariance sampling.

variance

float, optional

Variance of the model.

adaptive_subsampling

bool, default: True

Adaptive subsampling flag for full covariance calculation.

covariance_matrix

numpy.ndarray (pyrocko.guts_array.Array), optional

Cached covariance matrix, see covariance_matrix.

Model assumptions

kite.covariance.modelCovarianceExponential(distance, a, b)[source]

Exponential function model to approximate a positive-definite covariance

We assume the following simple covariance model to describe the empirical noise observations:

\[cov(d) = c \cdot e^{\frac{-d}{b}}\]
Parameters
  • distance (float or numpy.ndarray) – Distance between

  • a (float) – Linear model coefficient

  • b (float) – Exponential model coefficient

Returns

Covariance at distance

Return type

numpy.ndarray

kite.covariance.modelCovarianceExponentialCosine(distance, a, b, c, d)[source]

Exponential function model to approximate a positive-definite covariance

We assume the following simple covariance model to describe the empirical noise observations:

\[\begin{split}cov(d) = c \\cdot e^{\\frac{-d}{b}} \\cdot \cos{\\frac{d-c}{d}}\end{split}\]
Parameters
  • distance (float or numpy.ndarray) – Distance between

  • a (float) – Linear model coefficient

  • b (float) – Exponential model coefficient

  • c (float) – Cosinus distance correction

  • c – Cosinus coefficient

Returns

Covariance at distance

Return type

numpy.ndarray

kite.covariance.modelPowerspec(k, beta, D)[source]

Exponential linear model to estimate a log-linear power spectrum

We assume the following log-linear model for a measured power spectrum:

\[pow(k) = \frac{k^\beta}{D}\]
Parameters
  • k (float or numpy.ndarray) – Wavenumber

  • a (float) – Exponential model factor

  • b (float) – Fractional model factor