Generalized sampling is a numerically stable framework for obtaining reconstructions of signals in different bases and frames from their samples. For example, one can use wavelet bases for reconstruction given frequency measurements.

In this paper, we will introduce a carefully documented toolbox for performing generalized sampling in Julia. Julia is a new language for technical computing with focus on performance, which is ideally suited to handle the large size problems often encountered in generalized sampling. The toolbox provides specialized solutions for the setup of Fourier bases and wavelets.

The performance of the toolbox is compared to existing implementations of generalized sampling in MATLAB.

Generalized sampling [

The theory of generalized sampling does not restrict the type of bases to consider, but the applications have focused on Fourier bases and multiscale representations like wavelets.

Mathematically, samples of a function _{n}_{n}_{∈N} consists of inner products {〈_{n}_{n}_{∈N}. In generalized sampling we want to use these samples to estimate the inner products {〈_{n}_{n}_{∈N}, where {_{n}_{n}_{∈N} is another basis for _{n}_{n}_{∈N} is used for reconstructing

An especific example that will be used later is where 〈_{m}

that has been adapted to the interval

Other Daubechies scaling functions are at least continuous, but cannot be written in closed form in terms of elementary functions.

The actual computation of the reconstruction coefficients
_{j}, s_{i}

For large

We have two goals with the GeneralizedSampling package: It should be fast and easy to use. We have therefore put effort into providing only a few necessary high-level functions and hiding the lower level details. In essence, we need to compute a representation of a change of basis matrix

An example included in the package is the reconstruction of a truncated cosine (with inspiration from [

A truncated cosine and approximations with Haar and Daubechies 4 scaling functions using different scales. It may be difficult to see the graph of the original truncated cosine.

A reconstruction like this can be computed in only a few lines of code once the frequency measurements are available. For the truncated cosine we have a closed-form expression for the Fourier transform:

In Julia this can be implemented in the following way with a ternary operator (the line break is only for typographical purposes):

With the Fourier transform defined the reconstruction can be computed as follows:

To compute the reconstruction with the Daubechies 4 scaling functions, one simply has to replace “

As an example of reconstruction of 2D functions/images, consider reconstruction of a simulated brain made with the Matlab [^{2} frequency measurements we reconstruct 256^{2} Daubechies 4 scaling functions and the result is seen in Fig.

From 512 × 512 frequency measurements (not shown) 256 × 256 Daubechies 4 scaling functions are computed and used to produce this image. The image is evaluated at scale 10, i.e., in 1024 × 1024 points.

Julia is a dynamic language with a fast JIT compiler. A consequence is that a function is compiled the first time it is called, causing an overhead in terms of time and memory. All subsequent calls are, however, without this compiling overhead. The runtimes reported in this section are

The examples were carried out on a normal laptop (2.60 GHz Intel Core i7, 8 GB RAM) running GNU/Linux and are summarized in Table

Runtime comparisons with the Matlab implementation from [

init | sol | Iter | ||||||
---|---|---|---|---|---|---|---|---|

Problem | Size | Program | time | mem | time | mem | ||

Uniform 1D | 8192 × 4096 | Matlab | 15.0 | 901 | 0.13 | 37 | 9 | 0.01 |

Julia | 0.34 | 6.8 | 0.04 | 0.7 | 12 | 0.003 | ||

Jitter 1D | 5463 × 2048 | Matlab | 10.0 | 627 | 0.28 | 59 | 20 | 0.13 |

Julia | 0.23 | 4.0 | 0.08 | 0.4 | 20 | 0.004 | ||

Uniform 2D | 512^{2} × 256^{2} |
Matlab | 0.96 | 59 | 5.2 | 1507 | 9 | 0.58 |

Julia | 0.15 | 146 | 17.6 | 17 | 16 | 1.10 | ||

Jitter 2D | 26 244 × 32^{2} |
Matlab | 104.8 | 6377 | 8.3 | 1270 | 50 | 0.17 |

Julia | 2.4 | 165 | 2.4 | 1.3 | 18 | 0.13 | ||

Spiral | 27 681 × 32^{2} |
Matlab | 107.1 | 6670 | 3.7 | 479 | 17 | 0.21 |

Julia | 2.8 | 261 | 2.2 | 1.4 | 16 | 0.14 |

In both packages the solution step is based on a conjugate gradient-like algorithm, where the computational cost is dominated by multiplication with the change of basis matrix and its adjoint. In Matlab the built-in

In GeneralizedSampling we have relied on Julia’s ability to write functions that modify their arguments in-place to drastically reduce the memory consumption in an iterative algorithm like conjugate gradients. Julia’s

Especially for fast runtimes it is not accurate to rely on timing a single run of a function (using e.g.

The total memory allocation of the solver is much smaller in the Julia code than in the Matlab code. In most cases, the memory allocation of the initialization code is smaller the Julia code than in the Matlab code, although this is probably not as important since it is only performed once.

Note that the “sol” times and number of iterations are not directly comparable, since the stopping criteria for the least squares solver may be different and [^{2}([0, 1]^{d}

An extensive automated test suite has been developed along with the main code and is tested with Travis CI and AppVeyor, allowing tests to be performed on both Unix-like platforms and Windows. All public functions have doc-strings and a detailed documentation is included with the package in reStructuredText that can be exported to a multitude of formats. An HTML version of the documentation is hosted on Read the Docs. All details and links are available on the Generalized-Sampling package’s GitHub page.

Furthermore, all examples have been tested by all authors for bugs and readability.

Any system capable of running Julia version 0.4 or higher.

The GeneralizedSampling package has been tested with Julia version 0.4 and 0.5.

None.

The GeneralizedSampling package depends on the following Julia packages:

NFFT (version 0.1.4)

Wavelets (version 0.5.1)

ArrayViews (version 0.6.4)

VoronoiCells (version 0.1.3)

IntervalWavelets (version 0.0.4)

Jacobsen, Robert Dahl (Aalborg University) – Development.

Nielsen, Morten (Aalborg University) – Testing and theoretical discussions.

Rasmussen, Morten Grud (Aalborg University) – Testing and theoretical discussions.

English.

The theory of generalized sampling is very generic with possible applications in many different fields, e.g. where the aquisition of samples are fixed by physical restrictions. The GeneralizedSampling package is easy to use and and does not require users to purchase additional software.

Some parts of the software are easily extended: If one wishes to change representation between a pair other than the Fourier and wavelet bases, this can easily be achieved by adding appropriate functions for computing representations of the appropriate change of basis matrix and for performing multiplication with this matrix. It is possible to contribute to the package through the GitHub interface.

The authors have no competing interests to declare.