Function invDeconv

src/api/c/deconvolution.cpp:266–305

template<typename InputType, typename RealType = float>
af_array invDeconv(const af_array in, const af_array ker, const float gamma,
                   const af_inverse_deconv_algo algo) {
    using T    = RealType;
    using CT   = typename std::conditional<std::is_same<T, double>::value,
                                         cdouble, cfloat>::type;
    auto input = castArray<T>(in);
    auto psf   = castArray<T>(ker);
    const dim4& idims = input.dims();
    const dim4& fdims = psf.dims();
    dim_t nElems      = 1;

    dim4 inUPad, psfUPad, inLPad, psfLPad, odims(1);

    auto index = calcPadInfo(inLPad, psfLPad, inUPad, psfUPad, odims, nElems,
                             idims, fdims);
    auto paddedIn =
        padArrayBorders<T>(input, inLPad, inUPad, AF_PAD_CLAMP_TO_EDGE);
    auto paddedPsf = padArrayBorders<T>(psf, psfLPad, psfUPad, AF_PAD_ZERO);
    const array<int, 4> shiftDims = {-int(fdims[0] / 2), -int(fdims[1] / 2), 0,
                                     0};

    auto shiftedPsf = shift(paddedPsf, shiftDims.data());

    auto I      = fft_r2c<CT, T>(paddedIn, BASE_DIM);
    auto P      = fft_r2c<CT, T>(shiftedPsf, BASE_DIM);
    auto Pc     = conj(P);
    auto numer  = arithOp<CT, af_mul_t>(I, Pc, I.dims());
    auto denom  = denominator(I, P, gamma, algo);
    auto absVal = detail::abs<T, CT>(denom);
    auto THRESH = createValueArray(I.dims(), scalar<T>(gamma));
    auto cond   = logicOp<T, af_ge_t>(absVal, THRESH, absVal.dims());
    auto val    = arithOp<CT, af_div_t>(numer, denom, numer.dims());

    select_scalar<CT, false>(val, cond, val, scalar<CT>(0.0));

    auto ival =
        fft_c2r<CT, T>(val, 1 / static_cast<double>(nElems), odims, BASE_DIM);

    return getHandle(createSubArray<T>(ival, index));
}

af_err af_inverse_deconv(af_array* out, const af_array in, const af_array psf,
                         const float gamma, const af_inverse_deconv_algo algo) {