## Topology-based Smoothing of 1D Functions

Reconstruct1D is a tool for smoothing one-dimensional data based on topological constraints. Consider some input data containing a number of local minima and maxima. Given a selected subset of these minima/maxima, Reconstruct1D will produce a smoothed version of the input data with the following properties:

- the output data contains only the selected subset of minima/maxima,
- the output data is smooth in the sense of C1 or C2 continuity (user parameter),
- the output data is as close as desired to the original data (user parameter).

When working with noisy input data, the tool is best paired with Persistence1D: a tool for extracting and filtering extrema of one-dimensional data. The tool can just as well be used to create smooth one-dimensional data with topological constraints from scratch as can be seen in an example below.

The package supports Python and Matlab.

The Matlab version of Reconstruct1D has been written by Yeara Kozlov and Tino Weinkauf, Max Planck Institute for Informatics, Saarbrücken, Germany. The Python version has been written by Tino Weinkauf, KTH Royal Institute of Technology, Stockholm, Sweden. You may use the code as you wish, it is in the public domain. If you find it useful, it would be nice to hear from you. Just drop us a line.

## Overview

The algorithm is available in Python and Matlab:

**Python:**This is a re-implementation of the older Matlab code. It is cleaner and more documented. Only a single import statement is required. The code depends on the Python packages numpy and qpsolvers. The latter provides a unified API to a larger number of Quadratic Programming solvers (e.g., quadprog, cvxpy, ecos, and others) of which at least one has to be available. Reconstruct1D assumes quadprog to be available by default. These dependencies can be quickly installed using pip:- pip install numpy
- pip install qpsolvers
- pip install quadprog

**Matlab:**the code can be used from within Matlab straightforwardly. The use of MOSEK can be enabled to achieve fast computation times.

Both version should be used in conjunction with their respective Python/Matlab Persistence1D counterpart.

### Input

- One-dimensional vector of float values. It is assumed that this represents the function values of a one-dimensional function.
- Indices of minima and maxima points to remain in the smoothed data.
- Degree of smoothing, either 'biharmonic' for C1 continuous output, or 'triharmonic' for C2 continuous output.
- Data weight as a float value between 0 and 1 determining how close the smooth data shall be to the input data.

### Output

- One-dimensional vector of float values representing the smoothed data. The vector has as many entries as the input vector.

## Example for Python

A simple function call suffices to get the smoothed data with as many sample points as were given in the original data.
All code is located in the **python** folder.

```
import matplotlib.pyplot as plt
import numpy as np
from persistence1d import RunPersistence
from reconstruct1d import RunReconstruction
#~ Generate data using sine function and different frequencies.
x = np.arange(600.0)
SineLowFreq = np.sin(x * 0.01 * np.pi)
SineMedFreq = 0.25 * np.sin(x * 0.01 * np.pi * 4.9)
SineHighFreq = 0.15 * np.sin(x * 0.01 * np.pi * 12.1)
InputData = SineLowFreq + SineMedFreq + SineHighFreq
#~ Compute the extrema of the given data and their persistence.
ExtremaAndPersistence = RunPersistence(InputData)
#~ Keep only those extrema with a persistence larger than 0.5.
FilteredIndices = [t[0] for t in ExtremaAndPersistence if t[1] > 0.5]
#~ This simple call is all you need to reconstruct a smooth function containing only the filtered extrema
SmoothData = RunReconstruction(InputData, FilteredIndices, 'biharmonic', 0.0000001)
#~ Plot original and smoothed data
fig, ax = plt.subplots()
ax.plot(range(0, len(InputData)), InputData, label="Original Data")
ax.plot(range(0, len(SmoothData)), SmoothData, label="Smooth Data")
ExtremaIndices = [t[0] for t in ExtremaAndPersistence]
ax.plot(ExtremaIndices, InputData[ExtremaIndices], marker='.', linestyle='')
ax.plot(FilteredIndices, InputData[FilteredIndices], marker='*', linestyle='')
ax.set(xlabel='data index', ylabel='data value')
ax.set_aspect(1.0/ax.get_data_ratio()*0.2)
ax.grid()
plt.legend()
plt.show()
```

Running the code from above yields this output:

## More Examples

The code for the following examples can be found in the file **examples_reconstruct1d.py** in the **python** folder.

## Example for Matlab

The matlab interface is just as convenient and the results are easily plotted.
The result of this script is one of the plots from above.
The Matlab scripts for smoothing the function are located in the **reconstruct1d** folder.

```
% Add Reconstruct1D folder to Matlab's path
addpath('..');
% Setup Persistence1D and MOSEK
setup_persistence1d();
turn_on_mosek();
% Load the data set
load '..\datasets\test_data.mat';
% Choose smoothness for the reconstructed function.
% 'biharmonic' smoothness guarantees that the reconstructed function is C1 smooth
% 'triharmonic' smoothness guarantees that the reconstructed function is C2 smooth
smoothness = 'biharmonic';
% Choose a threshold for persistence features
threshold = 0.2;
% The data term weight affects how closely the reconstructed function
% adheres to the data.
data_weight = 0.0000001;
x = reconstruct1d(data, threshold, smoothness, data_weight);
plot_reconstructed_data(data, x, smoothness, threshold, data_weight);
turn_off_mosek();
```

## Documentation

The download package comes with extensive documentation and examples.

## Download

The Matlab version of Reconstruct1D has been written by Yeara Kozlov and Tino Weinkauf, Max Planck Institute for Informatics, Saarbrücken, Germany. The Python version has been written by Tino Weinkauf, KTH Royal Institute of Technology, Stockholm, Sweden. If you find it useful, it would be nice to hear from you. Just drop us a line.