Transfer Function

Plot

Regression

Base

follows Gary’s TTF.m in iris_mt_scratch/egbert_codes-20210121T193218Z-001/egbert_codes/matlabPrototype_10-13-20/TF/classes

2021-07-02 Removed some prototype methods intended to edit specific station data when many stations are being processed. MMT methods to be addressed later.

2021-07-20: Addressing Issue #12. If we are going to use xarray it is tempting to use the frequency band centers as the axis for the arrays here, rather than simple integer indexing. This has the advantage of making the data structures more explicit and self describing. We can also continue to use integer indices to assign and access tf values if needed. However, one concern I have is that if we use floating point numbers for the frequencies (or periods) we run the risk of machine roundoff error giving problems down stream. One way around this is to add a .band_centers() method to FrequencyBands() which will provide is a list of band centers and then rather than access by the band center, we can use an access method that gets us the frequencies between the band_egdes, which will be a unique frequency … however, if we use overlapping bands this will in general get complicated. However, for an MT TF representation, we do not in general use overlapping bands. A reasonably general, and simple solution is to make FrequencyBands support an iterable of bands, accessable by integer position, and by center frequency.

class aurora.transfer_function.base.TransferFunction(tf_header, frequency_bands, **kwargs)

Bases: object

Class to contain transfer function array. 2021-07-21: adding processing_config so that it lives with this object. This will facilitate the writing of z-files if we use a band setup file that is not of EMTF style

Parameters: TF : numpy array

array of transfer functions: TF(Nout, Nin, Nperiods)

Tnumpy array

list of periods

Headertransfer_function_header.TransferFunctionHeader() object.

TF header contains local site header, remote site header if appropriate, and information about estimation approach???

cov_ss_invnumpy array

inverse signal power matrix. aka Cov_SS in EMTF matlab codes

cov_nnnumpy array

noise covariance matrix: aka Cov_NN in EMTF matlab codes

num_segmentsinteger array?

Number of samples used to estimate TF for each band, and for each output channel (might be different for different channels)

R2xarray.DataArray

multiple coherence for each output channel / band

FullCovboolean

true if full covariance is provided

properties (Dependent) StdErr % standard errors of TF components, same size and order as TF NBands freqs % inverse of period Nout Nin

Attributes

maximum_period

minimum_period

num_bands

temporary function to allow access to old property num_bands used in

num_channels_in

num_channels_out

periods

tf

Methods

set_tf(regression_estimator, period)	This sets TF elements for one band, using contents of TRegression object.
standard_error()	TODO: make this a property that returns self._standard_error so it doesn't compute every time you call it.

frequency_index
from_emtf_zfile
period_index
tf_to_df
to_emtf_zfile

frequency_index(frequency)

from_emtf_zfile()

property maximum_period

property minimum_period

property num_bands

temporary function to allow access to old property num_bands used in the matlab codes for initialization Returns num_bands : int

a count of the frequency bands associated with the TF

property num_channels_in

property num_channels_out

period_index(period)

property periods

set_tf(regression_estimator, period)

This sets TF elements for one band, using contents of TRegression object. This version assumes there are estimates for Nout output channels

standard_error()

TODO: make this a property that returns self._standard_error so it doesn’t compute every time you call it. Returns ——-

property tf

tf_to_df()

to_emtf_zfile()

EMTF Z File Helpers

These methods can possibly be moved under mt_metadata, or mth5

They extract info needed to setup emtf_z files.

aurora.transfer_function.emtf_z_file_helpers.clip_bands_from_z_file(z_path, n_bands_clip, output_z_path=None, n_sensors=5)

This function takes a z_file and clips periods off the end of it. It can come in handy sometimes – specifically for manipulating matlab results of synthetic data.

Parameters

z_path: Path or str

path to the z_file to read in and clip periods from

n_periods_clip: integer

how many periods to clip from the end of the zfile

overwrite: bool

whether to overwrite the zfile or rename it

n_sensors

aurora.transfer_function.emtf_z_file_helpers.get_default_orientation_block(n_ch=5)

Helper function used when working with matlab structs which do not have enough info to make headers Returns ——-

aurora.transfer_function.emtf_z_file_helpers.make_orientation_block_of_z_file(run_obj, channel_list=None)

Replicates emtz z-file metadata about orientation like this: 1 0.00 0.00 tes Hx 2 90.00 0.00 tes Hy 3 0.00 0.00 tes Hz 4 0.00 0.00 tes Ex 5 90.00 0.00 tes Ey

based on this fortran snippet:

write(3, 115) k, orient(1, k), orient(2, k), stname(1: 3), chid(k)

format(i5, 1x, f8.2, 1x, f8.2, 1x, a3, 2x, a6) #Fortran Format Parameters ———- run_obj

aurora.transfer_function.emtf_z_file_helpers.merge_tf_collection_to_match_z_file(aux_data, tf_collection)

Currently this is only used for the synthtetic test, but maybe useful for other tests. Given data from a z_file, and a tf_collection, the tf_collection may have several TF estimates at the same frequency from multiple decimation levels. This tries to make a single array as a function of period for all rho and phi Parameters ———- aux_data tf_collection

Iter Control

follows Gary’s IterControl.m in iris_mt_scratch/egbert_codes-20210121T193218Z-001/egbert_codes/matlabPrototype_10-13-20/TF/classes

class aurora.transfer_function.iter_control.IterControl(max_number_of_iterations=10, max_number_of_redescending_iterations=2, **kwargs)

Bases: object

Attributes

continue_redescending

correction_factor

TODO: This is an RME specific property.

Methods

converged(b, b0)

Parameters

property continue_redescending

converged(b, b0)

Parameters

bcomplex-valued numpy array

the most recent regression estimate

b0complex-valued numpy array

The previous regression estimate

Returns

converged: bool

True of the regression has terminated, False otherwise

Notes:

The variable maximum_change finds the maximum amplitude component of the vector

1-b/b0. Looking at the formula, one might want to cast this instead as

1 - abs(b/b0), however, that will be insensitive to phase changes in b,

which is complex valued. The way it is coded np.max(np.abs(1 - b / b0)) is

correct as it stands.

property correction_factor

TODO: This is an RME specific property. Suggest move r0, u0 and this method into an RME-config class.

See notes on usage in transfer_function.regression.helper_functions.rme_beta

Returns

correction_factorfloat

correction factor used for scaling sigma_squared

Transfer Function Collection

Placeholder class. Will probably evolve structure in time. This is a container to hold: 1. TransferFunctionHeader 2. Dictionary of TransferFunction Objects

Note that a single transfer function object is associated with a station, which we call the “local_station”. In a database of TFs we could add a column for local_station and one for reference station.

class aurora.transfer_function.transfer_function_collection.TransferFunctionCollection(**kwargs)

Bases: object

Attributes

channel_list

local_station_id

TODO: make this take the station_id directly from the header

number_of_decimation_levels

total_number_of_channels

total_number_of_frequencies

Methods

relabel_merged_decimation_levels_for_export()	This method was specifcally related to issue #93, but may not be needed afterall Returns -------
rho_phi_plot([show, aux_data, xy_or_yx, ...])	One-off plotting method intended only for the synthetic test data for aurora dev Parameters ---------- show aux_data xy_or_yx ttl_str
write_emtf_z_file(z_file_path[, run_obj, ...])	Could probably move this into EMTFUtils() class Based on EMTF/T/wrt_z.f

check_all_channels_present
get_merged_dict

property channel_list

check_all_channels_present()

get_merged_dict()

property local_station_id

TODO: make this take the station_id directly from the header Returns ——-

property number_of_decimation_levels

relabel_merged_decimation_levels_for_export()

This method was specifcally related to issue #93, but may not be needed afterall Returns ——-

rho_phi_plot(show=True, aux_data=None, xy_or_yx='xy', ttl_str='', x_axis_fontsize=25, y_axis_fontsize=25, ttl_fontsize=16, markersize=10, rho_ylims=[10, 1000], phi_ylims=[0, 90], **kwargs)

One-off plotting method intended only for the synthetic test data for aurora dev Parameters ———- show aux_data xy_or_yx ttl_str

property total_number_of_channels

property total_number_of_frequencies

write_emtf_z_file(z_file_path, run_obj=None, orientation_strs=None)

Could probably move this into EMTFUtils() class Based on EMTF/T/wrt_z.f

Issues to review: This seems to insist that channels be ordered: Hx, Hy, Hz, Ex, Ey

z_file_path : Path or str

Sample output for a band: period : 4.65455 decimation level 1 freq. band from 25 to 30 number of data point 2489 sampling freq. 1.000 Hz

Transfer Functions

0.2498E+00 0.1966E-03 0.3859E-04 0.2519E+00

-0.1458E-01 -0.2989E-01 -0.7283E+01 -0.7313E+01

0.7311E+01 0.7338E+01 -0.4087E-01 -0.1031E-01

Inverse Coherent Signal Power Matrix

0.3809E-07 -0.6261E-18

-0.3095E-09 0.4505E-09 0.3764E-07 0.7792E-17 Residual Covariance

0.3639E+02 0.0000E+00

-0.2604E+03 0.2280E+03 0.3090E+05 0.0000E+00

0.2483E+03 0.2688E+03 0.2660E+03 -0.6791E+03 0.3161E+05 0.0000E+00

Transfer Function Header

follows Gary’s TFHeader.m iris_mt_scratch/egbert_codes-20210121T193218Z-001/egbert_codes/matlabPrototype_10-13-20/TF/classes

class aurora.transfer_function.transfer_function_header.TransferFunctionHeader(**kwargs)

Bases: object

class for storing metadata for a TF estimate

This class should inherit the metadata from the remote and reference stations. As of 2021-07-20 the class functions with only the station_id and channel_id values.

<See Issue #41> </See Issue #41>

Attributes

local_channels

local_station

local_station_id

num_input_channels

num_output_channels

reference_station_id

property local_channels

property local_station

property local_station_id

property num_input_channels

property num_output_channels

property reference_station_id

TTFZ

follows Gary’s TTFZ.m in iris_mt_scratch/egbert_codes-20210121T193218Z-001/egbert_codes/matlabPrototype_10-13-20/TF/classes

class aurora.transfer_function.TTFZ.TTFZ(*args, **kwargs)

Bases: aurora.transfer_function.base.TransferFunction

subclass to support some more MT impedance specficic functions – initially just apparent resistivity and pbase for diagonal elements + rotation/fixed coordinate system

properties rho rho_se phi phi_se

Attributes

maximum_period

minimum_period

num_bands

temporary function to allow access to old property num_bands used in

num_channels_in

num_channels_out

periods

tf

Methods

apparent_resistivity([units])	ap_res(...) : computes app.
set_tf(regression_estimator, period)	This sets TF elements for one band, using contents of TRegression object.
standard_error()	TODO: make this a property that returns self._standard_error so it doesn't compute every time you call it.

frequency_index
from_emtf_zfile
period_index
tf_to_df
to_emtf_zfile

apparent_resistivity(units='SI')

ap_res(…) : computes app. res., phase, errors, given imped., cov. %USAGE: [rho,rho_se,ph,ph_se] = ap_res(z,sig_s,sig_e,periods) ; % Z = array of impedances (from Z_***** file) % sig_s = inverse signal covariance matrix (from Z_****** file) % sig_e = residual covariance matrix (from Z_****** file) % periods = array of periods (sec) Returns ——-