How to resolve the algorithm Deconvolution/1D step by step in the Python programming language

Published on 12 May 2024 09:40 PM
#Python

How to resolve the algorithm Deconvolution/1D step by step in the Python programming language

Table of Contents

Problem Statement

The convolution of two functions

F

{\displaystyle {\mathit {F}}}

and

H

{\displaystyle {\mathit {H}}}

of an integer variable is defined as the function

G

{\displaystyle {\mathit {G}}}

satisfying for all integers

n

{\displaystyle {\mathit {n}}}

. Assume

F ( n )

{\displaystyle F(n)}

can be non-zero only for

0

{\displaystyle 0}

≤

n

{\displaystyle {\mathit {n}}}

≤

|

F

|

{\displaystyle |{\mathit {F}}|}

, where

|

F

|

{\displaystyle |{\mathit {F}}|}

is the "length" of

F

{\displaystyle {\mathit {F}}}

, and similarly for

G

{\displaystyle {\mathit {G}}}

and

H

{\displaystyle {\mathit {H}}}

, so that the functions can be modeled as finite sequences by identifying

f

0

,

f

1

,

f

2

, …

{\displaystyle f_{0},f_{1},f_{2},\dots }

with

F ( 0 ) , F ( 1 ) , F ( 2 ) , …

{\displaystyle F(0),F(1),F(2),\dots }

, etc. Then for example, values of

|

F

|

= 6

{\displaystyle |{\mathit {F}}|=6}

and

|

H

|

= 5

{\displaystyle |{\mathit {H}}|=5}

would determine the following value of

g

{\displaystyle {\mathit {g}}}

by definition. We can write this in matrix form as: or For this task, implement a function (or method, procedure, subroutine, etc.) deconv to perform deconvolution (i.e., the inverse of convolution) by constructing and solving such a system of equations represented by the above matrix

A

{\displaystyle A}

for

h

{\displaystyle {\mathit {h}}}

given

f

{\displaystyle {\mathit {f}}}

and

g

{\displaystyle {\mathit {g}}}

.

h = [-8,-9,-3,-1,-6,7] f = [-3,-6,-1,8,-6,3,-1,-9,-9,3,-2,5,2,-2,-7,-1] g = [24,75,71,-34,3,22,-45,23,245,25,52,25,-67,-96,96,31,55,36,29,-43,-7]

Let's start with the solution:

Step by Step solution about How to resolve the algorithm Deconvolution/1D step by step in the Python programming language

Let's begin by breaking down the functions and their roles in the code:

ToReducedRowEchelonForm(M):

  • This function takes a matrix M as input and attempts to convert it into its reduced row echelon form (RREF).
  • RREF is a specific representation of a matrix where each row has exactly one leading 1 (the leftmost nonzero element in the row) and all other elements in the same column are zero.
  • The algorithm works by iteratively scanning the matrix's rows and columns, performing row operations (row swaps and scaling) to achieve the RREF.

pmtx(mtx):

  • This is a utility function that prints a matrix mtx in a readable format.

convolve(f, h):

  • This function performs a convolution operation between two sequences f and h.
  • Convolution is a mathematical operation that combines two sequences to produce a third sequence.
  • In this case, the convolution result is stored in a new sequence g.

deconvolve(g, f):

  • This function performs a deconvolution operation, which is the inverse of convolution.
  • Given a sequence g and a kernel f, it aims to find a sequence h such that g is the convolution of f with h.

In the main part of the code, we have:

  • Two sequences, h and f, with some predefined values.

  • A sequence g that is calculated by applying convolution on f and h.

  • The code then performs various tests:

  • It asserts that convolve(f, h) equals g, confirming that the convolution operation is working correctly.

  • It asserts that deconvolve(g, f) equals h, indicating that the deconvolution operation can recover the original sequence h.

  • Finally, it prints the results of applying convolution and deconvolution using the same sequences h, f, and g for different versions of the convolution and deconvolution operations:

  • The first version uses pure Python code with the convolve() and deconvolve() functions.

  • The second version uses NumPy's FFT-based convolution and deconvolution functions.

  • This demonstrates that the results obtained from both methods are the same, adding confidence in the correctness of the implemented algorithms.

Source code in the python programming language

def ToReducedRowEchelonForm( M ):
    if not M: return
    lead = 0
    rowCount = len(M)
    columnCount = len(M[0])
    for r in range(rowCount):
        if lead >= columnCount:
            return
        i = r
        while M[i][lead] == 0:
            i += 1
            if i == rowCount:
                i = r
                lead += 1
                if columnCount == lead:
                    return
        M[i],M[r] = M[r],M[i]
        lv = M[r][lead]
        M[r] = [ mrx / lv for mrx in M[r]]
        for i in range(rowCount):
            if i != r:
                lv = M[i][lead]
                M[i] = [ iv - lv*rv for rv,iv in zip(M[r],M[i])]
        lead += 1
    return M
 
def pmtx(mtx):
    print ('\n'.join(''.join(' %4s' % col for col in row) for row in mtx))
 
def convolve(f, h):
    g = [0] * (len(f) + len(h) - 1)
    for hindex, hval in enumerate(h):
        for findex, fval in enumerate(f):
            g[hindex + findex] += fval * hval
    return g

def deconvolve(g, f):
    lenh = len(g) - len(f) + 1
    mtx = [[0 for x in range(lenh+1)] for y in g]
    for hindex in range(lenh):
        for findex, fval in enumerate(f):
            gindex = hindex + findex
            mtx[gindex][hindex] = fval
    for gindex, gval in enumerate(g):        
        mtx[gindex][lenh] = gval
    ToReducedRowEchelonForm( mtx )
    return [mtx[i][lenh] for i in range(lenh)]  # h

if __name__ == '__main__':
    h = [-8,-9,-3,-1,-6,7]
    f = [-3,-6,-1,8,-6,3,-1,-9,-9,3,-2,5,2,-2,-7,-1]
    g = [24,75,71,-34,3,22,-45,23,245,25,52,25,-67,-96,96,31,55,36,29,-43,-7]
    assert convolve(f,h) == g
    assert deconvolve(g, f) == h


import numpy

h = [-8,-9,-3,-1,-6,7]
f = [-3,-6,-1,8,-6,3,-1,-9,-9,3,-2,5,2,-2,-7,-1]
g = [24,75,71,-34,3,22,-45,23,245,25,52,25,-67,-96,96,31,55,36,29,-43,-7]

# https://stackoverflow.com/questions/14267555/find-the-smallest-power-of-2-greater-than-n-in-python

def shift_bit_length(x):
    return 1<<(x-1).bit_length()

def conv(a, b):
    p = len(a)
    q = len(b)
    n = p + q - 1
    r = shift_bit_length(n)
    y = numpy.fft.ifft(numpy.fft.fft(a,r) * numpy.fft.fft(b,r),r)
    return numpy.trim_zeros(numpy.around(numpy.real(y),decimals=6))

def deconv(a, b):
    p = len(a)
    q = len(b)
    n = p - q + 1
    r = shift_bit_length(max(p, q))
    y = numpy.fft.ifft(numpy.fft.fft(a,r) / numpy.fft.fft(b,r), r)
    return numpy.trim_zeros(numpy.around(numpy.real(y),decimals=6))
    
# should return g
   
print(conv(h,f))

# should return h

print(deconv(g,f))

# should return f

print(deconv(g,h))

You may also check:How to resolve the algorithm Speech synthesis step by step in the FutureBasic programming language
You may also check:How to resolve the algorithm Best shuffle step by step in the Rust programming language
You may also check:How to resolve the algorithm Polymorphism step by step in the Java programming language
You may also check:How to resolve the algorithm Statistics/Basic step by step in the Go programming language
You may also check:How to resolve the algorithm Bitwise operations step by step in the FALSE programming language