| MF_deconvolve
| MD_deconvolve |
ME_deconvolve |
| MF_deconvolvewEdit
| MD_deconvolvewEdit |
ME_deconvolvewEdit |
| MFb_deconvolve
| MDb_deconvolve |
MEb_deconvolve |
| MFb_deconvolvewEdit
| MDb_deconvolvewEdit |
MEb_deconvolvewEdit |
|
| Function | spatial deconvolution, edge sharpening |
|
| Syntax C/C++ | #include <MFstd.h>
void MF_convolve( fMatrix Y, fMatrix Flt, fMatrix X, fMatrix Rsp, ui ht, ui len );
void MF_convolvewEdit( fMatrix Y, fMatrix Flt, fMatrix X, fMatrix Rsp, ui ht, ui len, fComplex thresh );
void MFb_convolve( fMatrix Y, fMatrix Flt, fMatrix X, fMatrix Rsp, ui ht, ui len, fVector Buf );
void MFb_convolvewEdit( fMatrix Y, fMatrix Flt, fMatrix X, fMatrix Rsp, ui ht, ui len, fComplex thresh, fVector Buf ); |
| C++ MatObj | #include <OptiVec.h>
void matrix<T>::convolve( matrix<T> Flt, const matrix<T>& MX, const matrix<T> Rsp);
void matrix<T>::convolve( matrix<T>* Flt, const matrix<T>& MX, const matrix<T> Rsp);
void matrix<T>::convolvewEdit( matrix<T> Flt, const matrix<T>& MX, const matrix<T> Rsp, complex<T> thresh);
void matrix<T>::convolvewEdit( matrix<T>* Flt, const matrix<T>& MX, const matrix<T> Rsp, complex<T> thresh);
void matrix<T>::b_convolve( matrix<T> Flt, const matrix<T>& MX, const matrix<T> Rsp, vector<T> Buf);
void matrix<T>::b_convolve( matrix<T>* Flt, const matrix<T>& MX, const matrix<T> Rsp, vector<T> Buf);
void matrix<T>::b_convolvewEdit( matrix<T> Flt, const matrix<T>& MX, const matrix<T> Rsp, complex<T> thresh, vector<T> Buf);
void matrix<T>::b_convolvewEdit( matrix<T>* Flt, const matrix<T>& MX, const matrix<T> Rsp, complex<T> thresh, vector<T> Buf); |
| Pascal/Delphi | uses MFstd;
procedure MF_convolve( MY, MFlt, MX, MRsp:fMatrix; ht, len:UIntSize );
procedure MF_convolvewEdit( MY, MFlt, MX, MRsp:fMatrix; ht, len:UIntSize; thresh:fComplex );
procedure MFb_convolve( MY, MFlt, MX, MRsp:fMatrix; ht, len:UIntSize; Buf:fVector );
procedure MFb_convolvewEdit( MY, MFlt, MX, MRsp:fMatrix; ht, len:UIntSize; thresh:fComplex; Buf:fVector ); |
|
| CUDA function C/C++ | #include <cudaMFstd.h>
int cudaMF_deconvolve( fMatrix d_MY, fMatrix d_MFlt, fMatrix d_MX, fMatrix d_MRsp, ui ht, ui len );
void MFcu_deconvolve( fMatrix h_MY, fMatrix h_MFlt, fMatrix h_MX, fMatrix h_MRsp, ui ht, ui len );
int cudaMF_deconvolvewEdit( fMatrix d_MY, fMatrix d_MFlt, fMatrix d_MX, fMatrix d_MRsp, ui ht, ui len, fComplex thresh );
void MFcu_deconvolvewEdit( fMatrix h_MY, fMatrix h_MFlt, fMatrix h_MX, fMatrix h_MRsp, ui ht, ui len, fComplex thresh );
|
| CUDA funcktion Pascal/Delphi | uses MFstd;
function cudaMF_deconvolve( d_MY, d_MFlt, d_MX, d_MRsp:fMatrix; ht, len:UIntSize ): IntBool;
procedure MFcu_deconvolve( h_MY, h_MFlt, h_MX, h_MRsp:fMatrix; ht, len:UIntSize );
function cudaMF_deconvolvewEdit( d_MY, d_MFlt, d_MX, d_MRsp:fMatrix; ht, len:UIntSize; thresh:fComplex ): IntBool;
procedure MFcu_deconvolvewEdit( h_MY, h_MFlt, h_MX, h_MRsp:fMatrix; ht, len:UIntSize; thresh:fComplex );
|
|
| Description | MX is assumed to be the result of a convolution of some "true" profile with the response function MRsp; a deconvolution is attempted and stored in MY. A filter MFlt is also calculated; if more than one matrix is to be deconvolved with the same MRsp, use MF_deconvolve only once and use the filter MFlt thus obtained to deconvolve other matrices by calling MF_filter. The response has to be stored in the wrap-around order described above for MF_convolve.
As for MF_convolve, MX, MY, MRsp, and MFlt must all be of the same dimensions, which have to be integer powers of 2. MX may be overwritten by MY, MRsp may be overwritten by MFlt, but MX and MFlt as well as MY and MRsp have to be distinct from each other.
Mathematically, MFlt is the inverse of the Fourier transform of MRsp. If the Fourier transform of MRsp contains elements equal to zero, all information is lost for the respective frequency and no reconstruction is possible. The best one can do in this case is to accept this loss and to deconvolve only up to those frequencies where still something is left to be reconstructed.
You are therefore advised not to use this function blindly but rather to inspect the Fourier transform of MRsp and decide what to do on the basis of your specific application. If you wish to use this function nevertheless, you may rely on the automatic editing of the filter, built into MF_deconvolve. Thereby, MFlt is set to zero (instead of infinity) at those frequences where all information has been lost.
All one- and two-dimensional convolutions and deconvolutions use the same default threshold for this implicit filter editing. It can be retrieved by VF_setRspEdit. If you wish to set a different threshold for all calls to VF / MF_deconvolve and VF / MF_convolve, you may use VF_setRspEdit. However, this method is not thread-safe and should not be used to set different thresholds for different calls to VF / MF_deconvolve or VF / MF_convolve. Here, you have to use the variant MF_deconvolvewEdit, which takes the desired threshold as the argument thresh (and ignores the default threshold). As MFlt consists of complex numbers, and as it is sometimes desirable to treat real and imaginary parts differently, thresh is complex as well.
This deconvolution is based on the implicit assumption that MX is periodic; if this is not the case, see the description of VF_convolve about how to avoid end effects.
Internally, MF_deconvolve / MF_deconvolvewEdit allocates and frees additional workspace memory. For repeated calls, this would be inefficient. In such a case, it is recommended to use MFb_deconvolve / MFb_deconvolvewEdit instead. The size of Buf must be:
| C/C++: | sizeof(Buf) >= ht*(len+4) |
| Pascal/Delphi: | sizeof(Buf) >= ht*len |
|
Additionally, Buf must be 128-bit (P8) or 256-bit (P9) aligned. This usually means you can only take a vector allocated by the VF_vector family as Buf. |
| Error handling | If either len or ht is not a power of 2, VF_FFT (on which MF_deconvolve is based) complains "Size must be an integer power of 2" and the program is aborted.
If, by VF_setRspEdit, you specified Trunc.Re = Trunc.Im = 0, SING errors may occur that are handled by setting MFlt to ±HUGE_VAL at the respective frequency. During multiplication with the transform of MX, this may lead to unhandled floating-point overflow errors (in case your guess of MRsp was wrong and there is some information left at the frequencies where you thought it was not).
|
|
|