This repository was adapted from the supplementary material of
Torres Mendonça, G. L., Pongratz, J., and Reick, C. H.: Identification of linear response functions from arbitrary perturbation experiments in the presence of noise – Part 1: Method development and toy model demonstration, Nonlin. Processes Geophys., 28, 501–532, https://doi.org/10.5194/npg-28-501-2021, 2021.
The original code used for the publication can be found at
http://hdl.handle.net/21.11116/0000-0008-0F02-6
The RFI algorithm identifies linear response functions using data from a perturbation and a control experiment. All needed input information is entered in RFI_user.py, which is also the run script that calls the main routine RFI_main.py. After setting up RFI_user.py, the identification of response functions is performed by simply executing
python RFI_user.py
In RFI_user.py, the USER INTERFACE allows for selecting different setups to recover response functions from a toy model as well as from input data provided in the directory input_files. Such setups are switched on and off by assigning 1 and 0 values to the variables "$var$_sel".
In the first section "Main" the variable mode_sel determines whether the algorithm will be applied to a single time series or to many time series of different lengths to estimate the linear regime in the response according to the technique introduced in subsection 2.2 of the Part 2 paper. In this section other important parameters are defined, such as the parameter that determines whether the data for the recovery is artificial data from the toy model or imported data, and the parameters m, tau_min and tau_max needed for the application of the algorithm. Note that the values for tau_min and tau_max are actually log10(tau_min) and log10(tau_max) because the time-scale domain is taken as equally distributed at a logarithmic scale.
If the response function should be recovered from imported data, the corresponding data information is entered in the section "Setting for importing data (if applicable)" along with the names and units of the variables of interest. Important in this section is also to define whether the control data (both from the forcing and from the response) have some trend. If so, by activating the respective trend_controlresp_sel and trend_controlforc_sel, the algorithm automatically extracts these trends. If these options are deactivated, the equilibrium value for the control simulation is taken as the mean value over the control time series. Since taking the log forcing f0 ln (f(t)/f0) instead of f(t) is in some cases advantageous (see Part 2 paper), the option log_forcing_sel allows for automatically transforming the forcing in this way.
If the response function should be recovered from toy model simulations, the respective parameters are defined in "Setting for artificial experiment (if applicable)". Apart from the parameters of the toy model, in this section also the level of noise and nonlinearity in the artificial data can be adjusted.
In "Setting for processing" some technical parameters of the RFI algorithm are defined. Here the most important parameters are monotonicity_sel, which defines whether the additional noise level adjustment (step 6 of the RFI algorithm in Fig. 1 of Part 1) is switched on; monotonicity_tol_sel, which defines whether a tolerance should be used when judging the monotonicity of the response function; i_max_plateau_sel, which defines that in cases where the plateau in the singular values spectrum is not found, the high-frequency noise is defined by the two last singular values that are usually small; lamb_zero_sel, which defines whether the algorithm will tolerate a recovery with lambda = 0 (which happens for instance when no solution to the Discrepancy-method equation is found); lamb_min and lamb_max, which define the minimum and maximum values (in a log10 scale) for the lambda vector used to iteratively solve the Discrepancy-method equation.
"Setting for saving output data" defines the name of the output directory within "output_files".
"Setting for plotting" defines a standard setting for the final plotting.
- input_files/ and output_files/ must be in the root directory
- input files should be in .txt format with data organized in columns separated by an empty space, where it is by default assumed that the first column contains the time array and the second column the array for the variable of interest (see input_file/input_example.txt)