In [43]:

    

import sympy
from sympy import init_printing
from sympy import MutableDenseMatrix
# init_printing(use_latex='png')

# Use IPython's Latex capabilities to avoid a bug in Sympy's matrix printing
from IPython.display import display, Math, Latex

import numpy as np
import pyamg
from scipy import sparse
from scipy import ndimage
from scipy.linalg import solve
import rasterio as rio
from geoblend.coefficient_matrix import create_coefficient_matrix

import matplotlib as mpl
import matplotlib.pyplot as plt

mpl.rcParams['image.interpolation'] = 'none'
mpl.rcParams['image.cmap'] = 'cubehelix'

%matplotlib inline

np.set_printoptions(linewidth=200)


    

In [44]:

    

display(Math("\int \int_{\Omega} | f - \mathbf{v} | ^{2} \quad where \quad f |_{\delta \Omega} = f^{*}_{\delta \Omega}"))


    

    






$$\int \int_{\Omega} | f - \mathbf{v} | ^{2} \quad where \quad f |_{\delta \Omega} = f^{*}_{\delta \Omega}$$

In [45]:

    

def display_matrix(m):
    """
    Helper function to display the coefficient
    matrix.
    """
    
    if type(m) == MutableDenseMatrix:
        tex = sympy.latex(m)
    elif type(m) == np.ndarray:
        tex = sympy.latex(sympy.Matrix(mask.astype(np.int8).tolist()))
    else:
        # Crazy-ness to get integer represention ...
        tex = sympy.latex(sympy.Matrix(m.toarray().astype(np.int8).tolist()))
    
    display(Math(tex))


    

In [46]:

    

def create_auxiliary_equation(mask):
    """
    Create a sympy representation of the auxiliary equation corresponding
    to the Poisson equation.
    
    This prevents me from tedious, error-proned, hand-written algebra.
    """
    
    height, width = mask.shape
    x_domain = range(width)
    y_domain = range(height)
    
    # Create symbolic variables representing the unknown image
    # and the known gradient that will be imposed on the unknown image.
    
    # Note: One variable per pixel. f, s, c are pixel-aligned.
    # f represents the unknown image
    # s represents the source image
    # c represents the boundary conditions
    f = [ sympy.symbols(' '.join([ 'f%d%d' % (i, j) for j in y_domain ])) for i in x_domain ]
    s = [ sympy.symbols(' '.join([ 's%d%d' % (i, j) for j in y_domain ])) for i in x_domain ]
    c = [ sympy.symbols(' '.join([ 'c%d%d' % (i, j) for j in y_domain ])) for i in x_domain ]
    
    # Create the variable that will comprise the auxiliary equation.
    equation = 0
    
    for j in y_domain:
        for i in x_domain:
            
            if mask[j][i] == 0:
                continue
            
            #
            # Define the index of the 4-connected neighbors
            #
            
            # For a pixel p with indices (i, j), the
            # 4-connected neighbors are defined by the
            # surrounding 3x3 pixel grid as:
            #
            # - N - 
            # W p E
            # - S - 
            
            nj, ni = j - 1, i
            sj, si = j + 1, i
            ej, ei = j, i + 1
            wj, wi = j, i - 1
            
            if nj in y_domain:
                
                if mask[nj][ni] == 0:
                    equation += (f[i][j] - c[ni][nj] - (s[i][j] - s[ni][nj])) ** 2
                else:
                    equation += (f[i][j] - f[ni][nj] - (s[i][j] - s[ni][nj])) ** 2
            
            if sj in y_domain:
                
                if mask[sj][si] == 0:
                    equation += (f[i][j] - c[si][sj] - (s[i][j] - s[si][sj])) ** 2
                else:
                    equation += (f[i][j] - f[si][sj] - (s[i][j] - s[si][sj])) ** 2

            if ei in x_domain:
                
                if mask[ej][ei] == 0:
                    equation += (f[i][j] - c[ei][ej] - (s[i][j] - s[ei][ej])) ** 2
                else:
                    equation += (f[i][j] - f[ei][ej] - (s[i][j] - s[ei][ej])) ** 2
            
            if wi in x_domain:
                
                if mask[wj][wi] == 0:
                    equation += (f[i][j] - c[wi][wj] - (s[i][j] - s[wi][wj])) ** 2
                else:
                    equation += (f[i][j] - f[wi][wj] - (s[i][j] - s[wi][wj])) ** 2
    
    return equation, f, s, c


    

In [47]:

    

mask1 = np.array([
  [0, 0, 0, 0, 0, 0, 0, 0, 0, 0],
  [0, 1, 1, 1, 1, 1, 1, 1, 1, 0],
  [0, 1, 1, 1, 1, 1, 1, 1, 1, 0],
  [0, 1, 1, 1, 1, 1, 1, 1, 1, 0],
  [0, 1, 1, 1, 1, 1, 1, 1, 1, 0],
  [0, 1, 1, 1, 1, 1, 1, 1, 1, 0],
  [0, 1, 1, 1, 1, 1, 1, 1, 1, 0],
  [0, 0, 0, 0, 0, 0, 0, 0, 0, 0]
], dtype=np.uint8)


    

In [140]:

    

mask2 = np.array([
  [0, 0, 0, 0, 0, 0, 0, 0, 0, 0],
  [0, 0, 1, 1, 0, 0, 0, 0, 0, 0],
  [0, 0, 1, 1, 1, 0, 0, 0, 0, 0],
  [0, 0, 1, 1, 1, 1, 1, 1, 0, 0],
  [0, 0, 0, 1, 1, 1, 1, 1, 1, 0],
  [0, 0, 0, 1, 1, 1, 1, 1, 0, 0],
  [0, 0, 0, 0, 1, 1, 0, 0, 0, 0],
  [0, 0, 0, 0, 0, 0, 0, 0, 0, 0]
], dtype=np.uint8)
mask2


    

    
Out[140]:





array([[0, 0, 0, 0, 0, 0, 0, 0, 0, 0],
       [0, 0, 1, 1, 0, 0, 0, 0, 0, 0],
       [0, 0, 1, 1, 1, 0, 0, 0, 0, 0],
       [0, 0, 1, 1, 1, 1, 1, 1, 0, 0],
       [0, 0, 0, 1, 1, 1, 1, 1, 1, 0],
       [0, 0, 0, 1, 1, 1, 1, 1, 0, 0],
       [0, 0, 0, 0, 1, 1, 0, 0, 0, 0],
       [0, 0, 0, 0, 0, 0, 0, 0, 0, 0]], dtype=uint8)

In [49]:

    

# with rio.drivers():
#     with rio.open("test2/20150705_014439_081f_analytic.tif") as src:
#         mask3 = src.read(4).astype(np.uint8)
#         mask3[np.nonzero(mask3)] = 1
# mask3


    

In [169]:

    

mask4 = np.array([
    [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0],
    [0, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0],
    [0, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0],
    [0, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0],
    [0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0],
    [0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0],
    [0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0],
    [0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0],
    [0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0],
    [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]
], dtype=np.uint8)


    

In [170]:

    

# eq1, f1, s1, c1 = create_auxiliary_equation(mask1)
# eq2, f2, s2, c2 = create_auxiliary_equation(mask2)
# eq3, f3, s3, c3 = create_auxiliary_equation(mask3)
eq4, f4, s4, c4 = create_auxiliary_equation(mask4)


    

In [171]:

    

def linearize_system(mask, equation, f):
    
    indices = np.nonzero(mask)
    return [
        sympy.diff(equation, f[idx[1]][idx[0]])
        for idx in zip(*indices)
    ]


    

In [172]:

    

# system1 = linearize_system(mask1, eq1, f1)
# system2 = linearize_system(mask2, eq2, f2)
# system3 = linearize_system(mask3, eq3, f3)
system4 = linearize_system(mask4, eq4, f4)


    

In [173]:

    

# mat1 = sympy.Matrix(system1)
# mat2 = sympy.Matrix(system2)
# mat3 = sympy.Matrix(system3)
mat4 = sympy.Matrix(system4)


    

In [174]:

    

# display_matrix(mat1)
# display_matrix(mat2)
# display_matrix(mat3)
# display_matrix(mat4)


    

In [175]:

    

def get_unknown_vector(mask, f):
    
    indices = np.nonzero(mask)
    return [
        f[idx[1]][idx[0]]
        for idx in zip(*indices)
    ]


    

In [177]:

    

# vec1 = get_unknown_vector(mask1, f1)
# vec2 = get_unknown_vector(mask2, f2)
# vec3 = get_unknown_vector(mask3, f3)
vec4 = get_unknown_vector(mask4, f4)


    

In [179]:

    

# coeff1 = mat1.jacobian(vec1)
# coeff2 = mat2.jacobian(vec2)
# coeff3 = mat3.jacobian(vec3)
coeff4 = mat4.jacobian(vec4)


    

In [181]:

    

# display_matrix(coeff1)
# display_matrix(coeff2)
# display_matrix(coeff3)
coeff4


    

    
Out[181]:





Matrix([
[12, -4,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[-4, 14, -4,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0, -4, 14, -4,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0, -4, 12,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[-4,  0,  0,  0, 14, -4,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0, -4,  0,  0, -4, 16, -4,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0, -4,  0,  0, -4, 16, -4,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0, -4,  0,  0, -4, 16, -4,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0, -4, 14, -4,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0, -4, 12,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0, -4,  0,  0,  0,  0,  0, 14, -4,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0, -4,  0,  0,  0,  0, -4, 16, -4,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0, -4, 16, -4,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0, -4, 16, -4,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0, -4, 16, -4,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0, -4, 16, -4,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4, 12,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0, 14, -4,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0, -4, 16, -4,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0, -4, 16, -4,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0, -4, 16, -4,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0, -4, 16, -4,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0, -4, 16, -4,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0, -4, 16, -4,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4, 14, -4,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4, 12,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0, 14, -4,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0, -4, 16, -4,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0, -4, 16, -4,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0, -4, 16, -4,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0, -4, 16, -4,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0, -4, 16, -4,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0, -4, 16, -4,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0, -4, 16, -4,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0, -4, 16, -4,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4, 12,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0, 14, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0, -4, 16, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0, -4, 16, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0, -4, 16, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0, -4, 16, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0, -4, 16, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0, -4, 16, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0, -4, 16, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0, -4, 16, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0, -4, 16, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4, 12,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, 14, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4, 16, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4, 16, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4, 16, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4, 16, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4, 16, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4, 16, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4, 16, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4, 16, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4, 16, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4, 16, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4, 12,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, 12, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4, 14, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4, 14, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4, 14, -4,  0,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4, 14, -4,  0,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4, 14, -4,  0,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4, 14, -4,  0,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4, 14, -4,  0,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4, 14, -4,  0,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4, 14, -4,  0,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4, 14, -4,  0],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4, 14, -4],
[ 0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, -4, 10]])

In [26]:

    

# a = np.array(coeff3.tolist())
# np.save('matrix.npy', a)


    

In [95]:

    

# mask3


    

In [150]:

    

def get_replacements(source, reference, mask, f, c, s):
    
    boundary_replacements = []
    field_replacements = []
    
    # Replace boundaries with values from the reference image
    
    height, width = mask.shape
    for j in range(height):
        for i in range(width):
            
            boundary_replacements.append((c[i][j], reference[j][i]))
            field_replacements.append((s[i][j], source[j][i]))

    return boundary_replacements + field_replacements


    

In [151]:

    

def get_null_replacements_for_unknowns(mask, f):
    """
    Zero out unknown values from the system to get
    a grasp on the RHS.
    """
    replacements = []
    
    height, width = mask.shape
    for j in range(height):
        for i in range(width):
            replacements.append((f[i][j], 0))

    return replacements


    

In [136]:

    

null_replacements = get_null_replacements_for_unknowns(mask1, f1)


    

In [139]:

    

print mask1
display_matrix(mat1.subs(null_replacements))


    

    




[[0 0 0 0 0 0 0 0 0 0]
 [0 1 1 1 1 1 1 1 1 0]
 [0 1 1 1 1 1 1 1 1 0]
 [0 1 1 1 1 1 1 1 1 0]
 [0 1 1 1 1 1 1 1 1 0]
 [0 1 1 1 1 1 1 1 1 0]
 [0 1 1 1 1 1 1 1 1 0]
 [0 0 0 0 0 0 0 0 0 0]]






    






$$\left[\begin{matrix}- 2 c_{01} - 2 c_{10} + 2 s_{01} + 2 s_{10} - 12 s_{11} + 4 s_{12} + 4 s_{21}\\- 2 c_{20} + 4 s_{11} + 2 s_{20} - 14 s_{21} + 4 s_{22} + 4 s_{31}\\- 2 c_{30} + 4 s_{21} + 2 s_{30} - 14 s_{31} + 4 s_{32} + 4 s_{41}\\- 2 c_{40} + 4 s_{31} + 2 s_{40} - 14 s_{41} + 4 s_{42} + 4 s_{51}\\- 2 c_{50} + 4 s_{41} + 2 s_{50} - 14 s_{51} + 4 s_{52} + 4 s_{61}\\- 2 c_{60} + 4 s_{51} + 2 s_{60} - 14 s_{61} + 4 s_{62} + 4 s_{71}\\- 2 c_{70} + 4 s_{61} + 2 s_{70} - 14 s_{71} + 4 s_{72} + 4 s_{81}\\- 2 c_{80} - 2 c_{91} + 4 s_{71} + 2 s_{80} - 12 s_{81} + 4 s_{82} + 2 s_{91}\\- 2 c_{02} + 2 s_{02} + 4 s_{11} - 14 s_{12} + 4 s_{13} + 4 s_{22}\\4 s_{12} + 4 s_{21} - 16 s_{22} + 4 s_{23} + 4 s_{32}\\4 s_{22} + 4 s_{31} - 16 s_{32} + 4 s_{33} + 4 s_{42}\\4 s_{32} + 4 s_{41} - 16 s_{42} + 4 s_{43} + 4 s_{52}\\4 s_{42} + 4 s_{51} - 16 s_{52} + 4 s_{53} + 4 s_{62}\\4 s_{52} + 4 s_{61} - 16 s_{62} + 4 s_{63} + 4 s_{72}\\4 s_{62} + 4 s_{71} - 16 s_{72} + 4 s_{73} + 4 s_{82}\\- 2 c_{92} + 4 s_{72} + 4 s_{81} - 14 s_{82} + 4 s_{83} + 2 s_{92}\\- 2 c_{03} + 2 s_{03} + 4 s_{12} - 14 s_{13} + 4 s_{14} + 4 s_{23}\\4 s_{13} + 4 s_{22} - 16 s_{23} + 4 s_{24} + 4 s_{33}\\4 s_{23} + 4 s_{32} - 16 s_{33} + 4 s_{34} + 4 s_{43}\\4 s_{33} + 4 s_{42} - 16 s_{43} + 4 s_{44} + 4 s_{53}\\4 s_{43} + 4 s_{52} - 16 s_{53} + 4 s_{54} + 4 s_{63}\\4 s_{53} + 4 s_{62} - 16 s_{63} + 4 s_{64} + 4 s_{73}\\4 s_{63} + 4 s_{72} - 16 s_{73} + 4 s_{74} + 4 s_{83}\\- 2 c_{93} + 4 s_{73} + 4 s_{82} - 14 s_{83} + 4 s_{84} + 2 s_{93}\\- 2 c_{04} + 2 s_{04} + 4 s_{13} - 14 s_{14} + 4 s_{15} + 4 s_{24}\\4 s_{14} + 4 s_{23} - 16 s_{24} + 4 s_{25} + 4 s_{34}\\4 s_{24} + 4 s_{33} - 16 s_{34} + 4 s_{35} + 4 s_{44}\\4 s_{34} + 4 s_{43} - 16 s_{44} + 4 s_{45} + 4 s_{54}\\4 s_{44} + 4 s_{53} - 16 s_{54} + 4 s_{55} + 4 s_{64}\\4 s_{54} + 4 s_{63} - 16 s_{64} + 4 s_{65} + 4 s_{74}\\4 s_{64} + 4 s_{73} - 16 s_{74} + 4 s_{75} + 4 s_{84}\\- 2 c_{94} + 4 s_{74} + 4 s_{83} - 14 s_{84} + 4 s_{85} + 2 s_{94}\\- 2 c_{05} + 2 s_{05} + 4 s_{14} - 14 s_{15} + 4 s_{16} + 4 s_{25}\\4 s_{15} + 4 s_{24} - 16 s_{25} + 4 s_{26} + 4 s_{35}\\4 s_{25} + 4 s_{34} - 16 s_{35} + 4 s_{36} + 4 s_{45}\\4 s_{35} + 4 s_{44} - 16 s_{45} + 4 s_{46} + 4 s_{55}\\4 s_{45} + 4 s_{54} - 16 s_{55} + 4 s_{56} + 4 s_{65}\\4 s_{55} + 4 s_{64} - 16 s_{65} + 4 s_{66} + 4 s_{75}\\4 s_{65} + 4 s_{74} - 16 s_{75} + 4 s_{76} + 4 s_{85}\\- 2 c_{95} + 4 s_{75} + 4 s_{84} - 14 s_{85} + 4 s_{86} + 2 s_{95}\\- 2 c_{06} - 2 c_{17} + 2 s_{06} + 4 s_{15} - 12 s_{16} + 2 s_{17} + 4 s_{26}\\- 2 c_{27} + 4 s_{16} + 4 s_{25} - 14 s_{26} + 2 s_{27} + 4 s_{36}\\- 2 c_{37} + 4 s_{26} + 4 s_{35} - 14 s_{36} + 2 s_{37} + 4 s_{46}\\- 2 c_{47} + 4 s_{36} + 4 s_{45} - 14 s_{46} + 2 s_{47} + 4 s_{56}\\- 2 c_{57} + 4 s_{46} + 4 s_{55} - 14 s_{56} + 2 s_{57} + 4 s_{66}\\- 2 c_{67} + 4 s_{56} + 4 s_{65} - 14 s_{66} + 2 s_{67} + 4 s_{76}\\- 2 c_{77} + 4 s_{66} + 4 s_{75} - 14 s_{76} + 2 s_{77} + 4 s_{86}\\- 2 c_{87} - 2 c_{96} + 4 s_{76} + 4 s_{85} - 12 s_{86} + 2 s_{87} + 2 s_{96}\end{matrix}\right]$$

In [85]:

    

source1 = np.array([
    [502, 527, 545, 517, 518, 492, 457, 562, 405, 420],
    [605, 512, 444, 473, 465, 496, 527, 445, 387, 397],
    [543, 446, 440, 393, 491, 472, 471, 417, 439, 371],
    [513, 476, 494, 448, 470, 491, 492, 443, 559, 514],
    [454, 487, 498, 471, 402, 484, 471, 377, 574, 452],
    [507, 478, 499, 484, 381, 372, 249, 333, 607, 410],
    [451, 497, 497, 392, 389, 476, 357, 366, 400, 464],
    [485, 517, 567, 531, 443, 324, 370, 408, 361, 464]
], dtype=np.uint16)

reference1 = np.array([
    [611, 573, 639, 564, 626, 588, 556, 503, 458, 461],
    [689, 559, 532, 550, 572, 601, 521, 466, 469, 437],
    [631, 530, 513, 504, 545, 516, 428, 444, 447, 430],
    [648, 566, 518, 514, 592, 537, 518, 468, 658, 559],
    [553, 587, 556, 544, 423, 574, 546, 452, 456, 387],
    [590, 598, 583, 564, 408, 389, 219, 498, 501, 479],
    [565, 572, 564, 436, 442, 638, 208, 382, 466, 455],
    [566, 545, 570, 507, 429, 378, 425, 474, 425, 466]
], dtype=np.uint16)


    

In [86]:

    

replacements1 = get_replacements(source1, reference1, mask1, f1, c1, s1)


    

In [87]:

    

replacements1


    

    
Out[87]:





[(c00, 611),
 (c10, 573),
 (c20, 639),
 (c30, 564),
 (c40, 626),
 (c50, 588),
 (c60, 556),
 (c70, 503),
 (c80, 458),
 (c90, 461),
 (c01, 689),
 (c11, 559),
 (c21, 532),
 (c31, 550),
 (c41, 572),
 (c51, 601),
 (c61, 521),
 (c71, 466),
 (c81, 469),
 (c91, 437),
 (c02, 631),
 (c12, 530),
 (c22, 513),
 (c32, 504),
 (c42, 545),
 (c52, 516),
 (c62, 428),
 (c72, 444),
 (c82, 447),
 (c92, 430),
 (c03, 648),
 (c13, 566),
 (c23, 518),
 (c33, 514),
 (c43, 592),
 (c53, 537),
 (c63, 518),
 (c73, 468),
 (c83, 658),
 (c93, 559),
 (c04, 553),
 (c14, 587),
 (c24, 556),
 (c34, 544),
 (c44, 423),
 (c54, 574),
 (c64, 546),
 (c74, 452),
 (c84, 456),
 (c94, 387),
 (c05, 590),
 (c15, 598),
 (c25, 583),
 (c35, 564),
 (c45, 408),
 (c55, 389),
 (c65, 219),
 (c75, 498),
 (c85, 501),
 (c95, 479),
 (c06, 565),
 (c16, 572),
 (c26, 564),
 (c36, 436),
 (c46, 442),
 (c56, 638),
 (c66, 208),
 (c76, 382),
 (c86, 466),
 (c96, 455),
 (c07, 566),
 (c17, 545),
 (c27, 570),
 (c37, 507),
 (c47, 429),
 (c57, 378),
 (c67, 425),
 (c77, 474),
 (c87, 425),
 (c97, 466),
 (s00, 502),
 (s10, 527),
 (s20, 545),
 (s30, 517),
 (s40, 518),
 (s50, 492),
 (s60, 457),
 (s70, 562),
 (s80, 405),
 (s90, 420),
 (s01, 605),
 (s11, 512),
 (s21, 444),
 (s31, 473),
 (s41, 465),
 (s51, 496),
 (s61, 527),
 (s71, 445),
 (s81, 387),
 (s91, 397),
 (s02, 543),
 (s12, 446),
 (s22, 440),
 (s32, 393),
 (s42, 491),
 (s52, 472),
 (s62, 471),
 (s72, 417),
 (s82, 439),
 (s92, 371),
 (s03, 513),
 (s13, 476),
 (s23, 494),
 (s33, 448),
 (s43, 470),
 (s53, 491),
 (s63, 492),
 (s73, 443),
 (s83, 559),
 (s93, 514),
 (s04, 454),
 (s14, 487),
 (s24, 498),
 (s34, 471),
 (s44, 402),
 (s54, 484),
 (s64, 471),
 (s74, 377),
 (s84, 574),
 (s94, 452),
 (s05, 507),
 (s15, 478),
 (s25, 499),
 (s35, 484),
 (s45, 381),
 (s55, 372),
 (s65, 249),
 (s75, 333),
 (s85, 607),
 (s95, 410),
 (s06, 451),
 (s16, 497),
 (s26, 497),
 (s36, 392),
 (s46, 389),
 (s56, 476),
 (s66, 357),
 (s76, 366),
 (s86, 400),
 (s96, 464),
 (s07, 485),
 (s17, 517),
 (s27, 567),
 (s37, 531),
 (s47, 443),
 (s57, 324),
 (s67, 370),
 (s77, 408),
 (s87, 361),
 (s97, 464)]

In [88]:

    

display_matrix(mat1.subs(replacements1))


    

    






$$\left[\begin{matrix}12 f_{11} - 4 f_{12} - 4 f_{21} - 2844.0\\- 4 f_{11} + 14 f_{21} - 4 f_{22} - 4 f_{31} - 704.0\\- 4 f_{21} + 14 f_{31} - 4 f_{32} - 4 f_{41} - 1508.0\\- 4 f_{31} + 14 f_{41} - 4 f_{42} - 4 f_{51} - 886.0\\- 4 f_{41} + 14 f_{51} - 4 f_{52} - 4 f_{61} - 1280.0\\- 4 f_{51} + 14 f_{61} - 4 f_{62} - 4 f_{71} - 1928.0\\- 4 f_{61} + 14 f_{71} - 4 f_{72} - 4 f_{81} - 788.0\\- 4 f_{71} + 12 f_{81} - 4 f_{82} - 1294.0\\- 4 f_{11} + 14 f_{12} - 4 f_{13} - 4 f_{22} - 708.0\\- 4 f_{12} - 4 f_{21} + 16 f_{22} - 4 f_{23} - 4 f_{32} + 68.0\\- 4 f_{22} - 4 f_{31} + 16 f_{32} - 4 f_{33} - 4 f_{42} + 1120.0\\- 4 f_{32} - 4 f_{41} + 16 f_{42} - 4 f_{43} - 4 f_{52} - 656.0\\- 4 f_{42} - 4 f_{51} + 16 f_{52} - 4 f_{53} - 4 f_{62} + 244.0\\- 4 f_{52} - 4 f_{61} + 16 f_{62} - 4 f_{63} - 4 f_{72} + 96.0\\- 4 f_{62} - 4 f_{71} + 16 f_{72} - 4 f_{73} - 4 f_{82} + 520.0\\- 4 f_{72} - 4 f_{81} + 14 f_{82} - 4 f_{83} - 812.0\\- 4 f_{12} + 14 f_{13} - 4 f_{14} - 4 f_{23} - 1226.0\\- 4 f_{13} - 4 f_{22} + 16 f_{23} - 4 f_{24} - 4 f_{33} - 456.0\\- 4 f_{23} - 4 f_{32} + 16 f_{33} - 4 f_{34} - 4 f_{43} + 144.0\\- 4 f_{33} - 4 f_{42} + 16 f_{43} - 4 f_{44} - 4 f_{53} - 192.0\\- 4 f_{43} - 4 f_{52} + 16 f_{53} - 4 f_{54} - 4 f_{63} - 184.0\\- 4 f_{53} - 4 f_{62} + 16 f_{63} - 4 f_{64} - 4 f_{73} - 368.0\\- 4 f_{63} - 4 f_{72} + 16 f_{73} - 4 f_{74} - 4 f_{83} + 292.0\\- 4 f_{73} - 4 f_{82} + 14 f_{83} - 4 f_{84} - 2092.0\\- 4 f_{13} + 14 f_{14} - 4 f_{15} - 4 f_{24} - 1208.0\\- 4 f_{14} - 4 f_{23} + 16 f_{24} - 4 f_{25} - 4 f_{34} - 164.0\\- 4 f_{24} - 4 f_{33} + 16 f_{34} - 4 f_{35} - 4 f_{44} - 208.0\\- 4 f_{34} - 4 f_{43} + 16 f_{44} - 4 f_{45} - 4 f_{54} + 792.0\\- 4 f_{44} - 4 f_{53} + 16 f_{54} - 4 f_{55} - 4 f_{64} - 800.0\\- 4 f_{54} - 4 f_{63} + 16 f_{64} - 4 f_{65} - 4 f_{74} - 1128.0\\- 4 f_{64} - 4 f_{73} + 16 f_{74} - 4 f_{75} - 4 f_{84} + 1252.0\\- 4 f_{74} - 4 f_{83} + 14 f_{84} - 4 f_{85} - 1734.0\\- 4 f_{14} + 14 f_{15} - 4 f_{16} - 4 f_{25} - 926.0\\- 4 f_{15} - 4 f_{24} + 16 f_{25} - 4 f_{26} - 4 f_{35} - 156.0\\- 4 f_{25} - 4 f_{34} + 16 f_{35} - 4 f_{36} - 4 f_{45} - 772.0\\- 4 f_{35} - 4 f_{44} + 16 f_{45} - 4 f_{46} - 4 f_{55} + 492.0\\- 4 f_{45} - 4 f_{54} + 16 f_{55} - 4 f_{56} - 4 f_{65} + 408.0\\- 4 f_{55} - 4 f_{64} + 16 f_{65} - 4 f_{66} - 4 f_{75} + 2148.0\\- 4 f_{65} - 4 f_{74} + 16 f_{75} - 4 f_{76} - 4 f_{85} + 1068.0\\- 4 f_{75} - 4 f_{84} + 14 f_{85} - 4 f_{86} - 3408.0\\- 4 f_{15} + 12 f_{16} - 4 f_{26} - 2348.0\\- 4 f_{16} - 4 f_{25} + 14 f_{26} - 4 f_{36} - 1412.0\\- 4 f_{26} - 4 f_{35} + 14 f_{36} - 4 f_{46} + 40.0\\- 4 f_{36} - 4 f_{45} + 14 f_{46} - 4 f_{56} - 422.0\\- 4 f_{46} - 4 f_{55} + 14 f_{56} - 4 f_{66} - 2300.0\\- 4 f_{56} - 4 f_{65} + 14 f_{66} - 4 f_{76} - 744.0\\- 4 f_{66} - 4 f_{75} + 14 f_{76} - 4 f_{86} - 896.0\\- 4 f_{76} - 4 f_{85} + 12 f_{86} - 1018.0\end{matrix}\right]$$

In [91]:

    

mat1.subs(replacements1)


    

    
Out[91]:





Matrix([
[                 12*f11 - 4*f12 - 4*f21 - 2844.0],
[         -4*f11 + 14*f21 - 4*f22 - 4*f31 - 704.0],
[        -4*f21 + 14*f31 - 4*f32 - 4*f41 - 1508.0],
[         -4*f31 + 14*f41 - 4*f42 - 4*f51 - 886.0],
[        -4*f41 + 14*f51 - 4*f52 - 4*f61 - 1280.0],
[        -4*f51 + 14*f61 - 4*f62 - 4*f71 - 1928.0],
[         -4*f61 + 14*f71 - 4*f72 - 4*f81 - 788.0],
[                -4*f71 + 12*f81 - 4*f82 - 1294.0],
[         -4*f11 + 14*f12 - 4*f13 - 4*f22 - 708.0],
[  -4*f12 - 4*f21 + 16*f22 - 4*f23 - 4*f32 + 68.0],
[-4*f22 - 4*f31 + 16*f32 - 4*f33 - 4*f42 + 1120.0],
[ -4*f32 - 4*f41 + 16*f42 - 4*f43 - 4*f52 - 656.0],
[ -4*f42 - 4*f51 + 16*f52 - 4*f53 - 4*f62 + 244.0],
[  -4*f52 - 4*f61 + 16*f62 - 4*f63 - 4*f72 + 96.0],
[ -4*f62 - 4*f71 + 16*f72 - 4*f73 - 4*f82 + 520.0],
[         -4*f72 - 4*f81 + 14*f82 - 4*f83 - 812.0],
[        -4*f12 + 14*f13 - 4*f14 - 4*f23 - 1226.0],
[ -4*f13 - 4*f22 + 16*f23 - 4*f24 - 4*f33 - 456.0],
[ -4*f23 - 4*f32 + 16*f33 - 4*f34 - 4*f43 + 144.0],
[ -4*f33 - 4*f42 + 16*f43 - 4*f44 - 4*f53 - 192.0],
[ -4*f43 - 4*f52 + 16*f53 - 4*f54 - 4*f63 - 184.0],
[ -4*f53 - 4*f62 + 16*f63 - 4*f64 - 4*f73 - 368.0],
[ -4*f63 - 4*f72 + 16*f73 - 4*f74 - 4*f83 + 292.0],
[        -4*f73 - 4*f82 + 14*f83 - 4*f84 - 2092.0],
[        -4*f13 + 14*f14 - 4*f15 - 4*f24 - 1208.0],
[ -4*f14 - 4*f23 + 16*f24 - 4*f25 - 4*f34 - 164.0],
[ -4*f24 - 4*f33 + 16*f34 - 4*f35 - 4*f44 - 208.0],
[ -4*f34 - 4*f43 + 16*f44 - 4*f45 - 4*f54 + 792.0],
[ -4*f44 - 4*f53 + 16*f54 - 4*f55 - 4*f64 - 800.0],
[-4*f54 - 4*f63 + 16*f64 - 4*f65 - 4*f74 - 1128.0],
[-4*f64 - 4*f73 + 16*f74 - 4*f75 - 4*f84 + 1252.0],
[        -4*f74 - 4*f83 + 14*f84 - 4*f85 - 1734.0],
[         -4*f14 + 14*f15 - 4*f16 - 4*f25 - 926.0],
[ -4*f15 - 4*f24 + 16*f25 - 4*f26 - 4*f35 - 156.0],
[ -4*f25 - 4*f34 + 16*f35 - 4*f36 - 4*f45 - 772.0],
[ -4*f35 - 4*f44 + 16*f45 - 4*f46 - 4*f55 + 492.0],
[ -4*f45 - 4*f54 + 16*f55 - 4*f56 - 4*f65 + 408.0],
[-4*f55 - 4*f64 + 16*f65 - 4*f66 - 4*f75 + 2148.0],
[-4*f65 - 4*f74 + 16*f75 - 4*f76 - 4*f85 + 1068.0],
[        -4*f75 - 4*f84 + 14*f85 - 4*f86 - 3408.0],
[                -4*f15 + 12*f16 - 4*f26 - 2348.0],
[        -4*f16 - 4*f25 + 14*f26 - 4*f36 - 1412.0],
[          -4*f26 - 4*f35 + 14*f36 - 4*f46 + 40.0],
[         -4*f36 - 4*f45 + 14*f46 - 4*f56 - 422.0],
[        -4*f46 - 4*f55 + 14*f56 - 4*f66 - 2300.0],
[         -4*f56 - 4*f65 + 14*f66 - 4*f76 - 744.0],
[         -4*f66 - 4*f75 + 14*f76 - 4*f86 - 896.0],
[                -4*f76 - 4*f85 + 12*f86 - 1018.0]])

In [154]:

    

# Numpy representation of the valued elements above
b1 = np.array([
    2844., 704., 1508., 886., 1280., 1928., 788., 1294.,
    708., -68., -1120., 656., -244., -96., -520., 812.,
    1226., 456., -144., 192.,   184.,   368.,  -292.,
    2092.,  1208.,   164.,   208.,  -792.,   800.,  1128.,
    -1252.,  1734.,   926.,   156.,   772.,  -492., -408.,
    -2148., -1068.,  3408.,  2348.,  1412.,   -40.,   422.,
    2300.,   744.,   896.,  1018.
])


    

In [155]:

    

A1 = np.array(coeff1.tolist(), dtype=np.int32)


    

In [156]:

    

blended1 = np.copy(reference1)
blended1[np.nonzero(mask1)] = solve(A1, b1)


    

In [166]:

    

blended1.tolist()


    

    
Out[166]:





[[611, 573, 639, 564, 626, 588, 556, 503, 458, 461],
 [689, 584, 517, 540, 537, 563, 583, 475, 426, 437],
 [631, 523, 511, 459, 554, 530, 521, 457, 480, 430],
 [648, 558, 564, 510, 527, 543, 538, 483, 595, 559],
 [553, 561, 561, 526, 452, 531, 514, 413, 596, 387],
 [590, 543, 553, 529, 423, 416, 292, 373, 645, 479],
 [565, 556, 538, 422, 419, 517, 401, 410, 436, 455],
 [566, 545, 570, 507, 429, 378, 425, 474, 425, 466]]

In [ ]:

    




    

In [157]:

    

f, ((ax1, ax2), (ax3, ax4)) = plt.subplots(2, 2, sharex=True, sharey=True, figsize=(16, 12))

ax1.set_title('Reference')
ax1.imshow(reference1)

ax2.set_title('Source')
ax2.imshow(source1)

ax3.set_title('Blended')
ax3.imshow(blended1)

ax4.set_title('Difference (source - blended)')
ax4.imshow(source1.astype(np.float) - blended1.astype(np.float))


    

    
Out[157]:





<matplotlib.image.AxesImage at 0x111b8e910>






    

In [158]:

    

# Data is the same as rectangular case. Mask is different above.

source2 = np.array([
    [502, 527, 545, 517, 518, 492, 457, 562, 405, 420],
    [605, 512, 444, 473, 465, 496, 527, 445, 387, 397],
    [543, 446, 440, 393, 491, 472, 471, 417, 439, 371],
    [513, 476, 494, 448, 470, 491, 492, 443, 559, 514],
    [454, 487, 498, 471, 402, 484, 471, 377, 574, 452],
    [507, 478, 499, 484, 381, 372, 249, 333, 607, 410],
    [451, 497, 497, 392, 389, 476, 357, 366, 400, 464],
    [485, 517, 567, 531, 443, 324, 370, 408, 361, 464]
], dtype=np.uint16)

reference2 = np.array([
    [611, 573, 639, 564, 626, 588, 556, 503, 458, 461],
    [689, 559, 532, 550, 572, 601, 521, 466, 469, 437],
    [631, 530, 513, 504, 545, 516, 428, 444, 447, 430],
    [648, 566, 518, 514, 592, 537, 518, 468, 658, 559],
    [553, 587, 556, 544, 423, 574, 546, 452, 456, 387],
    [590, 598, 583, 564, 408, 389, 219, 498, 501, 479],
    [565, 572, 564, 436, 442, 638, 208, 382, 466, 455],
    [566, 545, 570, 507, 429, 378, 425, 474, 425, 466]
], dtype=np.uint16)


    

In [159]:

    

replacements2 = get_replacements(source2, reference2, mask2, f2, c2, s2)


    

In [184]:

    

mat2.subs(replacements2)


    

    
Out[184]:





Matrix([
[                 12*f21 - 4*f22 - 4*f31 - 1958.0],
[                -4*f21 + 12*f31 - 4*f32 - 2636.0],
[        -4*f21 + 14*f22 - 4*f23 - 4*f32 - 1004.0],
[-4*f22 - 4*f31 + 16*f32 - 4*f33 - 4*f42 + 1120.0],
[                -4*f32 + 12*f42 - 4*f43 - 2742.0],
[                -4*f22 + 12*f23 - 4*f33 - 2672.0],
[ -4*f23 - 4*f32 + 16*f33 - 4*f34 - 4*f43 + 144.0],
[ -4*f33 - 4*f42 + 16*f43 - 4*f44 - 4*f53 - 192.0],
[        -4*f43 + 14*f53 - 4*f54 - 4*f63 - 1178.0],
[        -4*f53 + 14*f63 - 4*f64 - 4*f73 - 1182.0],
[                -4*f63 + 12*f73 - 4*f74 - 2092.0],
[        -4*f33 + 14*f34 - 4*f35 - 4*f44 - 1374.0],
[ -4*f34 - 4*f43 + 16*f44 - 4*f45 - 4*f54 + 792.0],
[ -4*f44 - 4*f53 + 16*f54 - 4*f55 - 4*f64 - 800.0],
[-4*f54 - 4*f63 + 16*f64 - 4*f65 - 4*f74 - 1128.0],
[-4*f64 - 4*f73 + 16*f74 - 4*f75 - 4*f84 + 1252.0],
[                        -4*f74 + 10*f84 - 4088.0],
[                -4*f34 + 12*f35 - 4*f45 - 2656.0],
[ -4*f35 - 4*f44 + 16*f45 - 4*f46 - 4*f55 + 492.0],
[ -4*f45 - 4*f54 + 16*f55 - 4*f56 - 4*f65 + 408.0],
[        -4*f55 - 4*f64 + 14*f65 - 4*f75 + 1516.0],
[                -4*f65 - 4*f74 + 12*f75 - 1312.0],
[                -4*f45 + 12*f46 - 4*f56 - 1300.0],
[                -4*f46 - 4*f55 + 12*f56 - 2478.0]])

In [161]:

    

b2 = np.array([
        1958.,  2636.,  1004., -1120.,  2742.,  2672.,  -144.,
        192.,  1178.,  1182.,  2092.,  1374.,  -792.,   800.,
        1128., -1252.,  4088.,  2656.,  -492.,  -408., -1516.,
        1312.,  1300.,  2478.])


    

In [162]:

    

A2 = np.array(coeff2.tolist(), dtype=np.int32)


    

In [167]:

    

blended2 = np.copy(reference2)
blended2[np.nonzero(mask2)] = solve(A2, b2)


    

In [168]:

    

blended2.tolist()


    

    
Out[168]:





[[611, 573, 639, 564, 626, 588, 556, 503, 458, 461],
 [689, 559, 513, 542, 572, 601, 521, 466, 469, 437],
 [631, 530, 508, 455, 550, 516, 428, 444, 447, 430],
 [648, 566, 558, 500, 509, 512, 495, 462, 658, 559],
 [553, 587, 556, 513, 429, 492, 463, 368, 556, 387],
 [590, 598, 583, 525, 397, 365, 215, 303, 501, 479],
 [565, 572, 564, 436, 394, 459, 208, 382, 466, 455],
 [566, 545, 570, 507, 429, 378, 425, 474, 425, 466]]

In [ ]:

    




    

In [164]:

    

f, ((ax1, ax2), (ax3, ax4)) = plt.subplots(2, 2, sharex=True, sharey=True, figsize=(16, 12))

ax1.set_title('Reference')
ax1.imshow(reference2)

ax2.set_title('Source')
ax2.imshow(source2)

ax3.set_title('Blended')
ax3.imshow(blended2)

ax4.set_title('Difference (source - blended)')
ax4.imshow(source2.astype(np.float) - blended2.astype(np.float))


    

    
Out[164]:





<matplotlib.image.AxesImage at 0x111fd8490>






    

In [182]:

    

source4 = np.array([
    [231, 199, 236, 235, 242, 243,   0,   0,   0,   0,   0,   0,   0,   0,   0,   0,   0,   0,   0,   0],
    [266, 223, 236, 244, 248, 243, 227,   0,   0,   0,   0,   0,   0,   0,   0,   0,   0,   0,   0,   0],
    [221, 195, 230, 243, 235, 220, 221, 221,   0,   0,   0,   0,   0,   0,   0,   0,   0,   0,   0,   0],
    [214, 200, 226, 228, 232, 228, 226, 217, 229,   0,   0,   0,   0,   0,   0,   0,   0,   0,   0,   0],
    [210, 199, 215, 213, 229, 231, 218, 208, 221, 235,   0,   0,   0,   0,   0,   0,   0,   0,   0,   0],
    [227, 220, 227, 219, 207, 222, 208, 224, 210, 215, 262,   0,   0,   0,   0,   0,   0,   0,   0,   0],
    [224, 215, 218, 212, 202, 194, 216, 207, 209, 249, 232, 252,   0,   0,   0,   0,   0,   0,   0,   0],
    [226, 227, 221, 201, 189, 206, 208, 225, 236, 222, 255, 244, 222,   0,   0,   0,   0,   0,   0,   0],
    [222, 217, 209, 202, 196, 206, 227, 222, 233, 226, 250, 217, 245, 240,   0,   0,   0,   0,   0,   0],
    [211, 212, 210, 203, 199, 220, 218, 218, 230, 258, 248, 234, 233, 251, 273,   0,   0,   0,   0,   0]
], dtype=np.uint16)

reference4 = np.array([
    [277, 277, 277, 270, 270, 270, 270, 270, 270, 262, 262, 262, 262, 262, 262, 286, 286, 286, 286, 286],
    [288, 288, 288, 274, 274, 274, 270, 270, 270, 262, 262, 262, 262, 262, 262, 286, 286, 286, 286, 286],
    [288, 288, 288, 274, 274, 274, 274, 274, 274, 251, 251, 251, 251, 251, 251, 281, 281, 281, 281, 281],
    [288, 288, 288, 274, 274, 274, 274, 274, 274, 251, 251, 251, 251, 251, 281, 281, 281, 281, 281, 281],
    [288, 288, 288, 274, 274, 274, 274, 274, 274, 251, 251, 251, 251, 251, 281, 281, 281, 281, 281, 281],
    [288, 288, 288, 274, 274, 274, 274, 274, 274, 251, 251, 251, 251, 251, 281, 281, 281, 281, 281, 281],
    [288, 288, 288, 274, 274, 274, 274, 274, 274, 251, 251, 251, 251, 251, 281, 281, 281, 281, 281, 281],
    [290, 290, 290, 271, 271, 271, 271, 271, 271, 272, 272, 272, 272, 272, 279, 279, 279, 279, 279, 279],
    [290, 290, 290, 271, 271, 271, 271, 271, 271, 272, 272, 272, 272, 272, 279, 279, 279, 279, 279, 279],
    [290, 290, 290, 271, 271, 271, 271, 271, 271, 272, 272, 272, 272, 272, 279, 279, 279, 279, 279, 279]
], dtype=np.uint16)


    

In [183]:

    

replacements4 = get_replacements(source4, reference4, mask4, f4, c4, s4)


    

In [185]:

    

mat4.subs(replacements4)


    

    
Out[185]:





Matrix([
[                     12*f11 - 4*f12 - 4*f21 - 1152.0],
[             -4*f11 + 14*f21 - 4*f22 - 4*f31 - 598.0],
[             -4*f21 + 14*f31 - 4*f32 - 4*f41 - 578.0],
[                    -4*f31 + 12*f41 - 4*f42 - 1178.0],
[             -4*f11 + 14*f12 - 4*f13 - 4*f22 - 252.0],
[      -4*f12 - 4*f21 + 16*f22 - 4*f23 - 4*f32 - 80.0],
[     -4*f22 - 4*f31 + 16*f32 - 4*f33 - 4*f42 - 140.0],
[      -4*f32 - 4*f41 + 16*f42 - 4*f43 - 4*f52 + 12.0],
[             -4*f42 + 14*f52 - 4*f53 - 4*f62 - 406.0],
[                    -4*f52 + 12*f62 - 4*f63 - 1060.0],
[             -4*f12 + 14*f13 - 4*f14 - 4*f23 - 468.0],
[     -4*f13 - 4*f22 + 16*f23 - 4*f24 - 4*f33 - 124.0],
[       -4*f23 - 4*f32 + 16*f33 - 4*f34 - 4*f43 + 8.0],
[      -4*f33 - 4*f42 + 16*f43 - 4*f44 - 4*f53 - 32.0],
[      -4*f43 - 4*f52 + 16*f53 - 4*f54 - 4*f63 - 12.0],
[      -4*f53 - 4*f62 + 16*f63 - 4*f64 - 4*f73 - 80.0],
[                    -4*f63 + 12*f73 - 4*f74 - 1064.0],
[             -4*f13 + 14*f14 - 4*f15 - 4*f24 - 402.0],
[      -4*f14 - 4*f23 + 16*f24 - 4*f25 - 4*f34 + 20.0],
[     -4*f24 - 4*f33 + 16*f34 - 4*f35 - 4*f44 + 156.0],
[     -4*f34 - 4*f43 + 16*f44 - 4*f45 - 4*f54 - 132.0],
[     -4*f44 - 4*f53 + 16*f54 - 4*f55 - 4*f64 - 108.0],
[       -4*f54 - 4*f63 + 16*f64 - 4*f65 - 4*f74 + 4.0],
[     -4*f64 - 4*f73 + 16*f74 - 4*f75 - 4*f84 + 192.0],
[             -4*f74 + 14*f84 - 4*f85 - 4*f94 - 572.0],
[                    -4*f84 + 12*f94 - 4*f95 - 2080.0],
[             -4*f14 + 14*f15 - 4*f16 - 4*f25 - 638.0],
[     -4*f15 - 4*f24 + 16*f25 - 4*f26 - 4*f35 - 144.0],
[      -4*f25 - 4*f34 + 16*f35 - 4*f36 - 4*f45 - 68.0],
[     -4*f35 - 4*f44 + 16*f45 - 4*f46 - 4*f55 + 176.0],
[     -4*f45 - 4*f54 + 16*f55 - 4*f56 - 4*f65 - 192.0],
[     -4*f55 - 4*f64 + 16*f65 - 4*f66 - 4*f75 + 192.0],
[     -4*f65 - 4*f74 + 16*f75 - 4*f76 - 4*f85 - 252.0],
[     -4*f75 - 4*f84 + 16*f85 - 4*f86 - 4*f95 + 116.0],
[    -4*f105 - 4*f85 - 4*f94 + 16*f95 - 4*f96 + 384.0],
[                   12*f105 - 4*f106 - 4*f95 - 2360.0],
[             -4*f15 + 14*f16 - 4*f17 - 4*f26 - 478.0],
[      -4*f16 - 4*f25 + 16*f26 - 4*f27 - 4*f36 + 12.0],
[      -4*f26 - 4*f35 + 16*f36 - 4*f37 - 4*f46 - 32.0],
[      -4*f36 - 4*f45 + 16*f46 - 4*f47 - 4*f56 - 24.0],
[     -4*f46 - 4*f55 + 16*f56 - 4*f57 - 4*f66 + 280.0],
[     -4*f56 - 4*f65 + 16*f66 - 4*f67 - 4*f76 - 188.0],
[     -4*f66 - 4*f75 + 16*f76 - 4*f77 - 4*f86 + 184.0],
[     -4*f76 - 4*f85 + 16*f86 - 4*f87 - 4*f96 + 264.0],
[    -4*f106 - 4*f86 - 4*f95 + 16*f96 - 4*f97 - 472.0],
[ -4*f105 + 16*f106 - 4*f107 - 4*f116 - 4*f96 + 360.0],
[                 -4*f106 + 12*f116 - 4*f117 - 2124.0],
[             -4*f16 + 14*f17 - 4*f18 - 4*f27 - 694.0],
[     -4*f17 - 4*f26 + 16*f27 - 4*f28 - 4*f37 - 116.0],
[      -4*f27 - 4*f36 + 16*f37 - 4*f38 - 4*f47 + 80.0],
[     -4*f37 - 4*f46 + 16*f47 - 4*f48 - 4*f57 + 196.0],
[     -4*f47 - 4*f56 + 16*f57 - 4*f58 - 4*f67 - 108.0],
[     -4*f57 - 4*f66 + 16*f67 - 4*f68 - 4*f77 + 168.0],
[     -4*f67 - 4*f76 + 16*f77 - 4*f78 - 4*f87 - 108.0],
[     -4*f77 - 4*f86 + 16*f87 - 4*f88 - 4*f97 - 220.0],
[    -4*f107 - 4*f87 - 4*f96 + 16*f97 - 4*f98 + 312.0],
[ -4*f106 + 16*f107 - 4*f108 - 4*f117 - 4*f97 - 288.0],
[-4*f107 - 4*f116 + 16*f117 - 4*f118 - 4*f127 - 120.0],
[                 -4*f117 + 12*f127 - 4*f128 - 1754.0],
[                    -4*f17 + 12*f18 - 4*f28 - 1152.0],
[             -4*f18 - 4*f27 + 14*f28 - 4*f38 - 526.0],
[             -4*f28 - 4*f37 + 14*f38 - 4*f48 - 540.0],
[             -4*f38 - 4*f47 + 14*f48 - 4*f58 - 500.0],
[             -4*f48 - 4*f57 + 14*f58 - 4*f68 - 470.0],
[             -4*f58 - 4*f67 + 14*f68 - 4*f78 - 740.0],
[             -4*f68 - 4*f77 + 14*f78 - 4*f88 - 474.0],
[             -4*f78 - 4*f87 + 14*f88 - 4*f98 - 608.0],
[            -4*f108 - 4*f88 - 4*f97 + 14*f98 - 372.0],
[          -4*f107 + 14*f108 - 4*f118 - 4*f98 - 756.0],
[         -4*f108 - 4*f117 + 14*f118 - 4*f128 - 158.0],
[         -4*f118 - 4*f127 + 14*f128 - 4*f138 - 792.0],
[                          -4*f128 + 10*f138 - 2564.0]])

In [186]:

    

b4 = np.array([1152., 598., 578., 1178., 252., 80., 140., -12., 406., 1060., 468., 124., -8., 32., 12., 80., 1064., 402., -20., -156., 132., 108., -4., -192., 572., 2080., 638., 144., 68., -176., 192., -192., 252., -116., -384., 2360., 478., -12., 32., 24., -280., 188., -184., -264., 472., -360., 2124., 694., 116., -80., -196., 108., -168., 108., 220., -312., 288., 120., 1754., 1152., 526., 540., 500., 470., 740., 474., 608., 372., 756., 158., 792., 2564.])


    

In [197]:

    

img4 = np.copy(reference4).astype(np.float64)
img4[np.nonzero(mask4)] = b4


    

In [198]:

    

img4


    

    
Out[198]:





array([[  277.,   277.,   277.,   270.,   270.,   270.,   270.,   270.,   270.,   262.,   262.,   262.,   262.,   262.,   262.,   286.,   286.,   286.,   286.,   286.],
       [  288.,  1152.,   598.,   578.,  1178.,   274.,   270.,   270.,   270.,   262.,   262.,   262.,   262.,   262.,   262.,   286.,   286.,   286.,   286.,   286.],
       [  288.,   252.,    80.,   140.,   -12.,   406.,  1060.,   274.,   274.,   251.,   251.,   251.,   251.,   251.,   251.,   281.,   281.,   281.,   281.,   281.],
       [  288.,   468.,   124.,    -8.,    32.,    12.,    80.,  1064.,   274.,   251.,   251.,   251.,   251.,   251.,   281.,   281.,   281.,   281.,   281.,   281.],
       [  288.,   402.,   -20.,  -156.,   132.,   108.,    -4.,  -192.,   572.,  2080.,   251.,   251.,   251.,   251.,   281.,   281.,   281.,   281.,   281.,   281.],
       [  288.,   638.,   144.,    68.,  -176.,   192.,  -192.,   252.,  -116.,  -384.,  2360.,   251.,   251.,   251.,   281.,   281.,   281.,   281.,   281.,   281.],
       [  288.,   478.,   -12.,    32.,    24.,  -280.,   188.,  -184.,  -264.,   472.,  -360.,  2124.,   251.,   251.,   281.,   281.,   281.,   281.,   281.,   281.],
       [  290.,   694.,   116.,   -80.,  -196.,   108.,  -168.,   108.,   220.,  -312.,   288.,   120.,  1754.,   272.,   279.,   279.,   279.,   279.,   279.,   279.],
       [  290.,  1152.,   526.,   540.,   500.,   470.,   740.,   474.,   608.,   372.,   756.,   158.,   792.,  2564.,   279.,   279.,   279.,   279.,   279.,   279.],
       [  290.,   290.,   290.,   271.,   271.,   271.,   271.,   271.,   271.,   272.,   272.,   272.,   272.,   272.,   279.,   279.,   279.,   279.,   279.,   279.]])

In [199]:

    

print img4.min(),img4.max()
print np.median(img4)


    

    




-384.0 2564.0
277.0

In [200]:

    

plt.imshow(img4, vmin=-400, vmax=3000)
plt.colorbar()


    

    
Out[200]:





<matplotlib.colorbar.Colorbar instance at 0x11400f368>






    

In [195]:

    

plt.plot(np.arange(0, b4.size), b4)


    

    
Out[195]:





[<matplotlib.lines.Line2D at 0x113646090>]






    

In [202]:

    

A4 = np.array(coeff4.tolist(), dtype=np.int32)


    

In [203]:

    

blended4 = np.copy(reference4)
blended4[np.nonzero(mask4)] = solve(A4, b4)


    

In [209]:

    

plt.imshow(blended4)
plt.colorbar()


    

    
Out[209]:





<matplotlib.colorbar.Colorbar instance at 0x1144070e0>






    

In [222]:

    

blended_analytic = np.copy(reference4)


    

In [224]:

    

blended_analytic[np.nonzero(mask4)] = np.array(coeff4.LUsolve(sympy.Matrix(b4.tolist())).tolist()).astype(np.float64).ravel()


    

In [225]:

    

plt.imshow(blended_analytic)
plt.colorbar()


    

    
Out[225]:





<matplotlib.colorbar.Colorbar instance at 0x1145d1050>






    

In [29]:

    

with rio.drivers():
    with rio.open('test2/landsat.tif') as src:
        reference3 = src.read(1)
    with rio.open('test2/20150705_014439_081f_analytic.tif') as src:
        source3 = src.read(1)
        field3 = (pyamg.gallery.poisson(source.shape) * source.ravel()).reshape(source.shape)


    

In [30]:

    

replacements3 = get_replacements(mask3, f3, g3, reference3, field3)


    

In [31]:

    

mat3.subs(replacements3)


    

    
Out[31]:





Matrix([
[                          12*f115 - 4*f116 - 4*f125 + 558.0],
[               -4*f115 + 14*f125 - 4*f126 - 4*f135 + 2240.0],
[               -4*f125 + 14*f135 - 4*f136 - 4*f145 - 1608.0],
[               -4*f135 + 14*f145 - 4*f146 - 4*f155 + 2982.0],
[               -4*f145 + 14*f155 - 4*f156 - 4*f165 - 7848.0],
[                  -4*f155 + 14*f165 - 4*f166 - 4*f175 - 2.0],
[               -4*f165 + 14*f175 - 4*f176 - 4*f185 + 2028.0],
[                        -4*f175 + 12*f185 - 4*f186 - 4646.0],
[                          -4*f106 + 12*f96 - 4*f97 - 4268.0],
[                 14*f106 - 4*f107 - 4*f116 - 4*f96 + 7272.0],
[      -4*f106 - 4*f115 + 16*f116 - 4*f117 - 4*f126 - 5292.0],
[      -4*f116 - 4*f125 + 16*f126 - 4*f127 - 4*f136 - 5580.0],
[      -4*f126 - 4*f135 + 16*f136 - 4*f137 - 4*f146 + 2128.0],
[       -4*f136 - 4*f145 + 16*f146 - 4*f147 - 4*f156 + 868.0],
[      -4*f146 - 4*f155 + 16*f156 - 4*f157 - 4*f166 - 1728.0],
[      -4*f156 - 4*f165 + 16*f166 - 4*f167 - 4*f176 + 5548.0],
[      -4*f166 - 4*f175 + 16*f176 - 4*f177 - 4*f186 - 2508.0],
[      -4*f176 - 4*f185 + 16*f186 - 4*f187 - 4*f196 + 2508.0],
[               -4*f186 + 14*f196 - 4*f197 - 4*f206 - 4612.0],
[                        -4*f196 + 12*f206 - 4*f207 - 3018.0],
[                            12*f87 - 4*f88 - 4*f97 - 2318.0],
[          -4*f107 - 4*f87 - 4*f96 + 16*f97 - 4*f98 + 3040.0],
[       -4*f106 + 16*f107 - 4*f108 - 4*f117 - 4*f97 - 5984.0],
[      -4*f107 - 4*f116 + 16*f117 - 4*f118 - 4*f127 + 4928.0],
[       -4*f117 - 4*f126 + 16*f127 - 4*f128 - 4*f137 + 636.0],
[      -4*f127 - 4*f136 + 16*f137 - 4*f138 - 4*f147 + 3172.0],
[      -4*f137 - 4*f146 + 16*f147 - 4*f148 - 4*f157 - 7108.0],
[      -4*f147 - 4*f156 + 16*f157 - 4*f158 - 4*f167 + 4396.0],
[      -4*f157 - 4*f166 + 16*f167 - 4*f168 - 4*f177 - 2000.0],
[      -4*f167 - 4*f176 + 16*f177 - 4*f178 - 4*f187 - 6932.0],
[      -4*f177 - 4*f186 + 16*f187 - 4*f188 - 4*f197 + 4988.0],
[      -4*f187 - 4*f196 + 16*f197 - 4*f198 - 4*f207 + 2172.0],
[      -4*f197 - 4*f206 + 16*f207 - 4*f208 - 4*f217 + 2556.0],
[                       -4*f207 + 12*f217 - 4*f218 - 15122.0],
[                            12*f78 - 4*f79 - 4*f88 + 2424.0],
[           -4*f78 - 4*f87 + 16*f88 - 4*f89 - 4*f98 - 3632.0],
[            -4*f108 - 4*f88 - 4*f97 + 16*f98 - 4*f99 + 80.0],
[        -4*f107 + 16*f108 - 4*f109 - 4*f118 - 4*f98 + 700.0],
[      -4*f108 - 4*f117 + 16*f118 - 4*f119 - 4*f128 - 2956.0],
[      -4*f118 - 4*f127 + 16*f128 - 4*f129 - 4*f138 + 3304.0],
[      -4*f128 - 4*f137 + 16*f138 - 4*f139 - 4*f148 - 9624.0],
[      -4*f138 - 4*f147 + 16*f148 - 4*f149 - 4*f158 + 8732.0],
[      -4*f148 - 4*f157 + 16*f158 - 4*f159 - 4*f168 + 3748.0],
[      -4*f158 - 4*f167 + 16*f168 - 4*f169 - 4*f178 - 5224.0],
[      -4*f168 - 4*f177 + 16*f178 - 4*f179 - 4*f188 + 5116.0],
[      -4*f178 - 4*f187 + 16*f188 - 4*f189 - 4*f198 - 3932.0],
[      -4*f188 - 4*f197 + 16*f198 - 4*f199 - 4*f208 + 2608.0],
[      -4*f198 - 4*f207 + 16*f208 - 4*f209 - 4*f218 - 1160.0],
[      -4*f208 - 4*f217 + 16*f218 - 4*f219 - 4*f228 + 7860.0],
[                       -4*f218 + 12*f228 - 4*f229 - 11346.0],
[                          -4*f610 + 12*f69 - 4*f79 + 1418.0],
[          -4*f69 - 4*f710 - 4*f78 + 16*f79 - 4*f89 - 8276.0],
[          -4*f79 - 4*f810 - 4*f88 + 16*f89 - 4*f99 + 4564.0],
[          -4*f109 - 4*f89 - 4*f910 - 4*f98 + 16*f99 - 536.0],
[       -4*f1010 - 4*f108 + 16*f109 - 4*f119 - 4*f99 + 644.0],
[     -4*f109 - 4*f1110 - 4*f118 + 16*f119 - 4*f129 + 3268.0],
[     -4*f119 - 4*f1210 - 4*f128 + 16*f129 - 4*f139 - 5244.0],
[    -4*f129 - 4*f1310 - 4*f138 + 16*f139 - 4*f149 + 10272.0],
[     -4*f139 - 4*f1410 - 4*f148 + 16*f149 - 4*f159 - 2264.0],
[    -4*f149 - 4*f1510 - 4*f158 + 16*f159 - 4*f169 - 13032.0],
[    -4*f159 - 4*f1610 - 4*f168 + 16*f169 - 4*f179 + 14436.0],
[     -4*f169 - 4*f1710 - 4*f178 + 16*f179 - 4*f189 - 6568.0],
[     -4*f179 - 4*f1810 - 4*f188 + 16*f189 - 4*f199 - 1848.0],
[     -4*f189 - 4*f1910 - 4*f198 + 16*f199 - 4*f209 + 2336.0],
[     -4*f199 - 4*f2010 - 4*f208 + 16*f209 - 4*f219 - 1676.0],
[      -4*f209 - 4*f2110 - 4*f218 + 16*f219 - 4*f229 + 484.0],
[              -4*f219 - 4*f2210 - 4*f228 + 14*f229 - 3444.0],
[                         12*f510 - 4*f511 - 4*f610 - 5770.0],
[       -4*f510 + 16*f610 - 4*f611 - 4*f69 - 4*f710 - 3252.0],
[       -4*f610 + 16*f710 - 4*f711 - 4*f79 - 4*f810 + 1492.0],
[        -4*f710 + 16*f810 - 4*f811 - 4*f89 - 4*f910 - 124.0],
[      -4*f1010 - 4*f810 + 16*f910 - 4*f911 - 4*f99 + 2400.0],
[    16*f1010 - 4*f1011 - 4*f109 - 4*f1110 - 4*f910 + 2184.0],
[  -4*f1010 + 16*f1110 - 4*f1111 - 4*f119 - 4*f1210 - 1512.0],
[  -4*f1110 + 16*f1210 - 4*f1211 - 4*f129 - 4*f1310 - 1672.0],
[  -4*f1210 + 16*f1310 - 4*f1311 - 4*f139 - 4*f1410 - 2456.0],
[  -4*f1310 + 16*f1410 - 4*f1411 - 4*f149 - 4*f1510 - 1540.0],
[  -4*f1410 + 16*f1510 - 4*f1511 - 4*f159 - 4*f1610 + 3992.0],
[  -4*f1510 + 16*f1610 - 4*f1611 - 4*f169 - 4*f1710 + 1328.0],
[  -4*f1610 + 16*f1710 - 4*f1711 - 4*f179 - 4*f1810 - 6176.0],
[  -4*f1710 + 16*f1810 - 4*f1811 - 4*f189 - 4*f1910 + 6176.0],
[  -4*f1810 + 16*f1910 - 4*f1911 - 4*f199 - 4*f2010 - 4096.0],
[  -4*f1910 + 16*f2010 - 4*f2011 - 4*f209 - 4*f2110 + 4392.0],
[  -4*f2010 + 16*f2110 - 4*f2111 - 4*f219 - 4*f2210 - 1096.0],
[  -4*f2110 + 16*f2210 - 4*f2211 - 4*f229 - 4*f2310 + 5200.0],
[                    -4*f2210 + 12*f2310 - 4*f2311 - 10474.0],
[               -4*f510 + 14*f511 - 4*f512 - 4*f611 + 4590.0],
[      -4*f511 - 4*f610 + 16*f611 - 4*f612 - 4*f711 + 5248.0],
[       -4*f611 - 4*f710 + 16*f711 - 4*f712 - 4*f811 + 348.0],
[       -4*f711 - 4*f810 + 16*f811 - 4*f812 - 4*f911 - 328.0],
[     -4*f1011 - 4*f811 - 4*f910 + 16*f911 - 4*f912 - 6784.0],
[  -4*f1010 + 16*f1011 - 4*f1012 - 4*f1111 - 4*f911 - 1164.0],
[  -4*f1011 - 4*f1110 + 16*f1111 - 4*f1112 - 4*f1211 - 952.0],
[ -4*f1111 - 4*f1210 + 16*f1211 - 4*f1212 - 4*f1311 + 4048.0],
[ -4*f1211 - 4*f1310 + 16*f1311 - 4*f1312 - 4*f1411 - 3052.0],
[ -4*f1311 - 4*f1410 + 16*f1411 - 4*f1412 - 4*f1511 + 3328.0],
[  -4*f1411 - 4*f1510 + 16*f1511 - 4*f1512 - 4*f1611 - 616.0],
[ -4*f1511 - 4*f1610 + 16*f1611 - 4*f1612 - 4*f1711 - 4768.0],
[ -4*f1611 - 4*f1710 + 16*f1711 - 4*f1712 - 4*f1811 + 3088.0],
[  -4*f1711 - 4*f1810 + 16*f1811 - 4*f1812 - 4*f1911 + 596.0],
[  -4*f1811 - 4*f1910 + 16*f1911 - 4*f1912 - 4*f2011 + 508.0],
[ -4*f1911 - 4*f2010 + 16*f2011 - 4*f2012 - 4*f2111 - 5088.0],
[ -4*f2011 - 4*f2110 + 16*f2111 - 4*f2112 - 4*f2211 + 2692.0],
[  -4*f2111 - 4*f2210 + 16*f2211 - 4*f2212 - 4*f2311 + 676.0],
[           -4*f2211 - 4*f2310 + 14*f2311 - 4*f2312 - 5018.0],
[              -4*f511 + 14*f512 - 4*f513 - 4*f612 - 12474.0],
[      -4*f512 - 4*f611 + 16*f612 - 4*f613 - 4*f712 + 5608.0],
[      -4*f612 - 4*f711 + 16*f712 - 4*f713 - 4*f812 - 5716.0],
[      -4*f712 - 4*f811 + 16*f812 - 4*f813 - 4*f912 + 2016.0],
[     -4*f1012 - 4*f812 - 4*f911 + 16*f912 - 4*f913 + 1424.0],
[   -4*f1011 + 16*f1012 - 4*f1013 - 4*f1112 - 4*f912 - 612.0],
[ -4*f1012 - 4*f1111 + 16*f1112 - 4*f1113 - 4*f1212 + 3112.0],
[ -4*f1112 - 4*f1211 + 16*f1212 - 4*f1213 - 4*f1312 - 1916.0],
[ -4*f1212 - 4*f1311 + 16*f1312 - 4*f1313 - 4*f1412 + 5216.0],
[ -4*f1312 - 4*f1411 + 16*f1412 - 4*f1413 - 4*f1512 - 4268.0],
[ -4*f1412 - 4*f1511 + 16*f1512 - 4*f1513 - 4*f1612 + 1832.0],
[  -4*f1512 - 4*f1611 + 16*f1612 - 4*f1613 - 4*f1712 + 944.0],
[ -4*f1612 - 4*f1711 + 16*f1712 - 4*f1713 - 4*f1812 + 1328.0],
[ -4*f1712 - 4*f1811 + 16*f1812 - 4*f1813 - 4*f1912 - 1452.0],
[ -4*f1812 - 4*f1911 + 16*f1912 - 4*f1913 - 4*f2012 + 5632.0],
[  -4*f1912 - 4*f2011 + 16*f2012 - 4*f2013 - 4*f2112 - 888.0],
[ -4*f2012 - 4*f2111 + 16*f2112 - 4*f2113 - 4*f2212 - 3892.0],
[ -4*f2112 - 4*f2211 + 16*f2212 - 4*f2213 - 4*f2312 + 2296.0],
[ -4*f2212 - 4*f2311 + 16*f2312 - 4*f2313 - 4*f2412 + 1784.0],
[                    -4*f2312 + 12*f2412 - 4*f2413 - 11922.0],
[                          12*f413 - 4*f414 - 4*f513 - 266.0],
[      -4*f413 - 4*f512 + 16*f513 - 4*f514 - 4*f613 + 2276.0],
[       -4*f513 - 4*f612 + 16*f613 - 4*f614 - 4*f713 - 584.0],
[        -4*f613 - 4*f712 + 16*f713 - 4*f714 - 4*f813 - 20.0],
[      -4*f713 - 4*f812 + 16*f813 - 4*f814 - 4*f913 + 3708.0],
[     -4*f1013 - 4*f813 - 4*f912 + 16*f913 - 4*f914 + 1580.0],
[  -4*f1012 + 16*f1013 - 4*f1014 - 4*f1113 - 4*f913 - 2848.0],
[  -4*f1013 - 4*f1112 + 16*f1113 - 4*f1114 - 4*f1213 + 852.0],
[ -4*f1113 - 4*f1212 + 16*f1213 - 4*f1214 - 4*f1313 - 2016.0],
[  -4*f1213 - 4*f1312 + 16*f1313 - 4*f1314 - 4*f1413 + 312.0],
[  -4*f1313 - 4*f1412 + 16*f1413 - 4*f1414 - 4*f1513 - 864.0],
[  -4*f1413 - 4*f1512 + 16*f1513 - 4*f1514 - 4*f1613 + 380.0],
[  -4*f1513 - 4*f1612 + 16*f1613 - 4*f1614 - 4*f1713 - 548.0],
[  -4*f1613 - 4*f1712 + 16*f1713 - 4*f1714 - 4*f1813 + 828.0],
[ -4*f1713 - 4*f1812 + 16*f1813 - 4*f1814 - 4*f1913 - 2000.0],
[ -4*f1813 - 4*f1912 + 16*f1913 - 4*f1914 - 4*f2013 - 6740.0],
[ -4*f1913 - 4*f2012 + 16*f2013 - 4*f2014 - 4*f2113 + 5160.0],
[ -4*f2013 - 4*f2112 + 16*f2113 - 4*f2114 - 4*f2213 + 3432.0],
[ -4*f2113 - 4*f2212 + 16*f2213 - 4*f2214 - 4*f2313 - 2104.0],
[ -4*f2213 - 4*f2312 + 16*f2313 - 4*f2314 - 4*f2413 + 2340.0],
[ -4*f2313 - 4*f2412 + 16*f2413 - 4*f2414 - 4*f2513 + 3332.0],
[                     -4*f2413 + 12*f2513 - 4*f2514 - 9422.0],
[               -4*f413 + 14*f414 - 4*f415 - 4*f514 - 8312.0],
[      -4*f414 - 4*f513 + 16*f514 - 4*f515 - 4*f614 + 7636.0],
[      -4*f514 - 4*f613 + 16*f614 - 4*f615 - 4*f714 - 5752.0],
[      -4*f614 - 4*f713 + 16*f714 - 4*f715 - 4*f814 + 3624.0],
[      -4*f714 - 4*f813 + 16*f814 - 4*f815 - 4*f914 - 4400.0],
[     -4*f1014 - 4*f814 - 4*f913 + 16*f914 - 4*f915 - 1648.0],
[  -4*f1013 + 16*f1014 - 4*f1015 - 4*f1114 - 4*f914 + 4804.0],
[ -4*f1014 - 4*f1113 + 16*f1114 - 4*f1115 - 4*f1214 - 1824.0],
[  -4*f1114 - 4*f1213 + 16*f1214 - 4*f1215 - 4*f1314 + 332.0],
[  -4*f1214 - 4*f1313 + 16*f1314 - 4*f1315 - 4*f1414 - 832.0],
[ -4*f1314 - 4*f1413 + 16*f1414 - 4*f1415 - 4*f1514 + 3876.0],
[ -4*f1414 - 4*f1513 + 16*f1514 - 4*f1515 - 4*f1614 - 3088.0],
[ -4*f1514 - 4*f1613 + 16*f1614 - 4*f1615 - 4*f1714 - 6744.0],
[ -4*f1614 - 4*f1713 + 16*f1714 - 4*f1715 - 4*f1814 + 5980.0],
[ -4*f1714 - 4*f1813 + 16*f1814 - 4*f1815 - 4*f1914 - 2760.0],
[ -4*f1814 - 4*f1913 + 16*f1914 - 4*f1915 - 4*f2014 + 3540.0],
[ -4*f1914 - 4*f2013 + 16*f2014 - 4*f2015 - 4*f2114 - 1124.0],
[  -4*f2014 - 4*f2113 + 16*f2114 - 4*f2115 - 4*f2214 - 832.0],
[  -4*f2114 - 4*f2213 + 16*f2214 - 4*f2215 - 4*f2314 - 152.0],
[  -4*f2214 - 4*f2313 + 16*f2314 - 4*f2315 - 4*f2414 - 768.0],
[ -4*f2314 - 4*f2413 + 16*f2414 - 4*f2415 - 4*f2514 + 2276.0],
[           -4*f2414 - 4*f2513 + 14*f2514 - 4*f2515 - 4468.0],
[                -4*f414 + 14*f415 - 4*f416 - 4*f515 - 752.0],
[      -4*f415 - 4*f514 + 16*f515 - 4*f516 - 4*f615 - 2508.0],
[      -4*f515 - 4*f614 + 16*f615 - 4*f616 - 4*f715 + 6352.0],
[      -4*f615 - 4*f714 + 16*f715 - 4*f716 - 4*f815 - 5236.0],
[      -4*f715 - 4*f814 + 16*f815 - 4*f816 - 4*f915 + 4364.0],
[     -4*f1015 - 4*f815 - 4*f914 + 16*f915 - 4*f916 - 1532.0],
[  -4*f1014 + 16*f1015 - 4*f1016 - 4*f1115 - 4*f915 - 3308.0],
[ -4*f1015 - 4*f1114 + 16*f1115 - 4*f1116 - 4*f1215 + 2196.0],
[  -4*f1115 - 4*f1214 + 16*f1215 - 4*f1216 - 4*f1315 + 132.0],
[ -4*f1215 - 4*f1314 + 16*f1315 - 4*f1316 - 4*f1415 - 4548.0],
[  -4*f1315 - 4*f1414 + 16*f1415 - 4*f1416 - 4*f1515 + 996.0],
[ -4*f1415 - 4*f1514 + 16*f1515 - 4*f1516 - 4*f1615 + 1640.0],
[ -4*f1515 - 4*f1614 + 16*f1615 - 4*f1616 - 4*f1715 + 8112.0],
[ -4*f1615 - 4*f1714 + 16*f1715 - 4*f1716 - 4*f1815 + 3580.0],
[-4*f1715 - 4*f1814 + 16*f1815 - 4*f1816 - 4*f1915 - 13720.0],
[-4*f1815 - 4*f1914 + 16*f1915 - 4*f1916 - 4*f2015 + 11800.0],
[-4*f1915 - 4*f2014 + 16*f2015 - 4*f2016 - 4*f2115 - 10448.0],
[ -4*f2015 - 4*f2114 + 16*f2115 - 4*f2116 - 4*f2215 + 1960.0],
[  -4*f2115 - 4*f2214 + 16*f2215 - 4*f2216 - 4*f2315 + 996.0],
[  -4*f2215 - 4*f2314 + 16*f2315 - 4*f2316 - 4*f2415 + 136.0],
[ -4*f2315 - 4*f2414 + 16*f2415 - 4*f2416 - 4*f2515 + 1540.0],
[           -4*f2415 - 4*f2514 + 14*f2515 - 4*f2516 + 1022.0],
[                       -4*f415 + 12*f416 - 4*f516 - 11540.0],
[     -4*f416 - 4*f515 + 16*f516 - 4*f517 - 4*f616 + 12460.0],
[      -4*f516 - 4*f615 + 16*f616 - 4*f617 - 4*f716 - 6468.0],
[      -4*f616 - 4*f715 + 16*f716 - 4*f717 - 4*f816 + 1420.0],
[      -4*f716 - 4*f815 + 16*f816 - 4*f817 - 4*f916 - 1320.0],
[      -4*f1016 - 4*f816 - 4*f915 + 16*f916 - 4*f917 + 940.0],
[  -4*f1015 + 16*f1016 - 4*f1017 - 4*f1116 - 4*f916 + 2856.0],
[ -4*f1016 - 4*f1115 + 16*f1116 - 4*f1117 - 4*f1216 - 1488.0],
[  -4*f1116 - 4*f1215 + 16*f1216 - 4*f1217 - 4*f1316 - 512.0],
[ -4*f1216 - 4*f1315 + 16*f1316 - 4*f1317 - 4*f1416 + 6816.0],
[ -4*f1316 - 4*f1415 + 16*f1416 - 4*f1417 - 4*f1516 + 2636.0],
[-4*f1416 - 4*f1515 + 16*f1516 - 4*f1517 - 4*f1616 - 10152.0],
[   -4*f1516 - 4*f1615 + 16*f1616 - 4*f1617 - 4*f1716 - 88.0],
[  -4*f1616 - 4*f1715 + 16*f1716 - 4*f1717 - 4*f1816 - 744.0],
[ -4*f1716 - 4*f1815 + 16*f1816 - 4*f1817 - 4*f1916 + 6008.0],
[ -4*f1816 - 4*f1915 + 16*f1916 - 4*f1917 - 4*f2016 - 5704.0],
[ -4*f1916 - 4*f2015 + 16*f2016 - 4*f2017 - 4*f2116 + 6252.0],
[  -4*f2016 - 4*f2115 + 16*f2116 - 4*f2117 - 4*f2216 + 256.0],
[ -4*f2116 - 4*f2215 + 16*f2216 - 4*f2217 - 4*f2316 - 1012.0],
[  -4*f2216 - 4*f2315 + 16*f2316 - 4*f2317 - 4*f2416 - 208.0],
[  -4*f2316 - 4*f2415 + 16*f2416 - 4*f2417 - 4*f2516 - 652.0],
[           -4*f2416 - 4*f2515 + 14*f2516 - 4*f2517 - 2494.0],
[                       -4*f516 + 12*f517 - 4*f617 - 13294.0],
[     -4*f517 - 4*f616 + 16*f617 - 4*f618 - 4*f717 + 10916.0],
[      -4*f617 - 4*f716 + 16*f717 - 4*f718 - 4*f817 + 1160.0],
[       -4*f717 - 4*f816 + 16*f817 - 4*f818 - 4*f917 - 692.0],
[     -4*f1017 - 4*f817 - 4*f916 + 16*f917 - 4*f918 - 1928.0],
[  -4*f1016 + 16*f1017 - 4*f1018 - 4*f1117 - 4*f917 + 2900.0],
[  -4*f1017 - 4*f1116 + 16*f1117 - 4*f1118 - 4*f1217 - 708.0],
[ -4*f1117 - 4*f1216 + 16*f1217 - 4*f1218 - 4*f1317 - 3064.0],
[ -4*f1217 - 4*f1316 + 16*f1317 - 4*f1318 - 4*f1417 - 4516.0],
[ -4*f1317 - 4*f1416 + 16*f1417 - 4*f1418 - 4*f1517 - 2916.0],
[ -4*f1417 - 4*f1516 + 16*f1517 - 4*f1518 - 4*f1617 + 9248.0],
[  -4*f1517 - 4*f1616 + 16*f1617 - 4*f1618 - 4*f1717 + 436.0],
[ -4*f1617 - 4*f1716 + 16*f1717 - 4*f1718 - 4*f1817 - 4804.0],
[ -4*f1717 - 4*f1816 + 16*f1817 - 4*f1818 - 4*f1917 + 2820.0],
[ -4*f1817 - 4*f1916 + 16*f1917 - 4*f1918 - 4*f2017 - 4004.0],
[  -4*f1917 - 4*f2016 + 16*f2017 - 4*f2018 - 4*f2117 - 676.0],
[ -4*f2017 - 4*f2116 + 16*f2117 - 4*f2118 - 4*f2217 + 2904.0],
[ -4*f2117 - 4*f2216 + 16*f2217 - 4*f2218 - 4*f2317 - 1112.0],
[  -4*f2217 - 4*f2316 + 16*f2317 - 4*f2318 - 4*f2417 - 344.0],
[  -4*f2317 - 4*f2416 + 16*f2417 - 4*f2418 - 4*f2517 - 768.0],
[                      -4*f2417 - 4*f2516 + 12*f2517 + 296.0],
[               -4*f617 + 14*f618 - 4*f619 - 4*f718 - 3646.0],
[        -4*f618 - 4*f717 + 16*f718 - 4*f719 - 4*f818 - 76.0],
[      -4*f718 - 4*f817 + 16*f818 - 4*f819 - 4*f918 - 2532.0],
[     -4*f1018 - 4*f818 - 4*f917 + 16*f918 - 4*f919 + 4052.0],
[  -4*f1017 + 16*f1018 - 4*f1019 - 4*f1118 - 4*f918 - 5504.0],
[ -4*f1018 - 4*f1117 + 16*f1118 - 4*f1119 - 4*f1218 + 4048.0],
[ -4*f1118 - 4*f1217 + 16*f1218 - 4*f1219 - 4*f1318 + 2656.0],
[  -4*f1218 - 4*f1317 + 16*f1318 - 4*f1319 - 4*f1418 - 576.0],
[ -4*f1318 - 4*f1417 + 16*f1418 - 4*f1419 - 4*f1518 + 1296.0],
[ -4*f1418 - 4*f1517 + 16*f1518 - 4*f1519 - 4*f1618 - 2692.0],
[ -4*f1518 - 4*f1617 + 16*f1618 - 4*f1619 - 4*f1718 - 5272.0],
[ -4*f1618 - 4*f1717 + 16*f1718 - 4*f1719 - 4*f1818 + 6176.0],
[ -4*f1718 - 4*f1817 + 16*f1818 - 4*f1819 - 4*f1918 + 2100.0],
[  -4*f1818 - 4*f1917 + 16*f1918 - 4*f1919 - 4*f2018 - 888.0],
[ -4*f1918 - 4*f2017 + 16*f2018 - 4*f2019 - 4*f2118 + 2768.0],
[ -4*f2018 - 4*f2117 + 16*f2118 - 4*f2119 - 4*f2218 - 5096.0],
[  -4*f2118 - 4*f2217 + 16*f2218 - 4*f2219 - 4*f2318 + 676.0],
[ -4*f2218 - 4*f2317 + 16*f2318 - 4*f2319 - 4*f2418 + 1960.0],
[           -4*f2318 - 4*f2417 + 14*f2418 - 4*f2419 - 1012.0],
[                        -4*f618 + 12*f619 - 4*f719 - 6810.0],
[      -4*f619 - 4*f718 + 16*f719 - 4*f720 - 4*f819 + 5216.0],
[      -4*f719 - 4*f818 + 16*f819 - 4*f820 - 4*f919 - 4208.0],
[     -4*f1019 - 4*f819 - 4*f918 + 16*f919 - 4*f920 + 1228.0],
[   -4*f1018 + 16*f1019 - 4*f1020 - 4*f1119 - 4*f919 - 944.0],
[ -4*f1019 - 4*f1118 + 16*f1119 - 4*f1120 - 4*f1219 - 1480.0],
[ -4*f1119 - 4*f1218 + 16*f1219 - 4*f1220 - 4*f1319 - 5428.0],
[-4*f1219 - 4*f1318 + 16*f1319 - 4*f1320 - 4*f1419 + 10132.0],
[ -4*f1319 - 4*f1418 + 16*f1419 - 4*f1420 - 4*f1519 - 5412.0],
[  -4*f1419 - 4*f1518 + 16*f1519 - 4*f1520 - 4*f1619 + 772.0],
[ -4*f1519 - 4*f1618 + 16*f1619 - 4*f1620 - 4*f1719 + 8408.0],
[ -4*f1619 - 4*f1718 + 16*f1719 - 4*f1720 - 4*f1819 - 7972.0],
[ -4*f1719 - 4*f1818 + 16*f1819 - 4*f1820 - 4*f1919 - 3684.0],
[ -4*f1819 - 4*f1918 + 16*f1919 - 4*f1920 - 4*f2019 + 7220.0],
[ -4*f1919 - 4*f2018 + 16*f2019 - 4*f2020 - 4*f2119 - 3128.0],
[ -4*f2019 - 4*f2118 + 16*f2119 - 4*f2120 - 4*f2219 - 2084.0],
[ -4*f2119 - 4*f2218 + 16*f2219 - 4*f2220 - 4*f2319 + 5384.0],
[ -4*f2219 - 4*f2318 + 16*f2319 - 4*f2320 - 4*f2419 - 4220.0],
[           -4*f2319 - 4*f2418 + 14*f2419 - 4*f2420 + 1576.0],
[               -4*f719 + 14*f720 - 4*f721 - 4*f820 + 2496.0],
[      -4*f720 - 4*f819 + 16*f820 - 4*f821 - 4*f920 + 6404.0],
[     -4*f1020 - 4*f820 - 4*f919 + 16*f920 - 4*f921 - 4876.0],
[  -4*f1019 + 16*f1020 - 4*f1021 - 4*f1120 - 4*f920 + 6292.0],
[ -4*f1020 - 4*f1119 + 16*f1120 - 4*f1121 - 4*f1220 - 1300.0],
[ -4*f1120 - 4*f1219 + 16*f1220 - 4*f1221 - 4*f1320 + 6480.0],
[-4*f1220 - 4*f1319 + 16*f1320 - 4*f1321 - 4*f1420 - 15804.0],
[-4*f1320 - 4*f1419 + 16*f1420 - 4*f1421 - 4*f1520 + 16064.0],
[ -4*f1420 - 4*f1519 + 16*f1520 - 4*f1521 - 4*f1620 - 9592.0],
[ -4*f1520 - 4*f1619 + 16*f1620 - 4*f1621 - 4*f1720 - 4708.0],
[ -4*f1620 - 4*f1719 + 16*f1720 - 4*f1721 - 4*f1820 + 4780.0],
[ -4*f1720 - 4*f1819 + 16*f1820 - 4*f1821 - 4*f1920 + 2660.0],
[ -4*f1820 - 4*f1919 + 16*f1920 - 4*f1921 - 4*f2020 - 9108.0],
[ -4*f1920 - 4*f2019 + 16*f2020 - 4*f2021 - 4*f2120 + 5032.0],
[  -4*f2020 - 4*f2119 + 16*f2120 - 4*f2121 - 4*f2220 + 896.0],
[ -4*f2120 - 4*f2219 + 16*f2220 - 4*f2221 - 4*f2320 - 2588.0],
[  -4*f2220 - 4*f2319 + 16*f2320 - 4*f2321 - 4*f2420 + 528.0],
[                     -4*f2320 - 4*f2419 + 12*f2420 - 2886.0],
[                       -4*f720 + 12*f721 - 4*f821 - 21098.0],
[     -4*f721 - 4*f820 + 16*f821 - 4*f822 - 4*f921 + 13004.0],
[     -4*f1021 - 4*f821 - 4*f920 + 16*f921 - 4*f922 - 6780.0],
[   -4*f1020 + 16*f1021 - 4*f1022 - 4*f1121 - 4*f921 + 156.0],
[ -4*f1021 - 4*f1120 + 16*f1121 - 4*f1122 - 4*f1221 - 3760.0],
[ -4*f1121 - 4*f1220 + 16*f1221 - 4*f1222 - 4*f1321 + 3120.0],
[ -4*f1221 - 4*f1320 + 16*f1321 - 4*f1322 - 4*f1421 + 2136.0],
[ -4*f1321 - 4*f1420 + 16*f1421 - 4*f1422 - 4*f1521 - 7184.0],
[-4*f1421 - 4*f1520 + 16*f1521 - 4*f1522 - 4*f1621 + 12224.0],
[ -4*f1521 - 4*f1620 + 16*f1621 - 4*f1622 - 4*f1721 - 2572.0],
[  -4*f1621 - 4*f1720 + 16*f1721 - 4*f1722 - 4*f1821 + 364.0],
[ -4*f1721 - 4*f1820 + 16*f1821 - 4*f1822 - 4*f1921 - 2100.0],
[ -4*f1821 - 4*f1920 + 16*f1921 - 4*f1922 - 4*f2021 + 7332.0],
[ -4*f1921 - 4*f2020 + 16*f2021 - 4*f2022 - 4*f2121 - 4172.0],
[ -4*f2021 - 4*f2120 + 16*f2121 - 4*f2122 - 4*f2221 + 1504.0],
[  -4*f2121 - 4*f2220 + 16*f2221 - 4*f2222 - 4*f2321 - 344.0],
[                     -4*f2221 - 4*f2320 + 12*f2321 - 1458.0],
[                        -4*f821 + 12*f822 - 4*f922 - 8810.0],
[    -4*f1022 - 4*f822 - 4*f921 + 16*f922 - 4*f923 + 13168.0],
[  -4*f1021 + 16*f1022 - 4*f1023 - 4*f1122 - 4*f922 - 2904.0],
[ -4*f1022 - 4*f1121 + 16*f1122 - 4*f1123 - 4*f1222 + 2640.0],
[ -4*f1122 - 4*f1221 + 16*f1222 - 4*f1223 - 4*f1322 - 4984.0],
[ -4*f1222 - 4*f1321 + 16*f1322 - 4*f1323 - 4*f1422 + 3144.0],
[ -4*f1322 - 4*f1421 + 16*f1422 - 4*f1423 - 4*f1522 - 1292.0],
[ -4*f1422 - 4*f1521 + 16*f1522 - 4*f1523 - 4*f1622 - 4876.0],
[ -4*f1522 - 4*f1621 + 16*f1622 - 4*f1623 - 4*f1722 + 3480.0],
[  -4*f1622 - 4*f1721 + 16*f1722 - 4*f1723 - 4*f1822 + 988.0],
[ -4*f1722 - 4*f1821 + 16*f1822 - 4*f1823 - 4*f1922 - 5808.0],
[ -4*f1822 - 4*f1921 + 16*f1922 - 4*f1923 - 4*f2022 + 1736.0],
[  -4*f1922 - 4*f2021 + 16*f2022 - 4*f2023 - 4*f2122 - 688.0],
[ -4*f2022 - 4*f2121 + 16*f2122 - 4*f2123 - 4*f2222 + 2088.0],
[                     -4*f2122 - 4*f2221 + 12*f2222 - 3568.0],
[              -4*f1023 - 4*f922 + 14*f923 - 4*f924 - 9810.0],
[  -4*f1022 + 16*f1023 - 4*f1024 - 4*f1123 - 4*f923 + 9048.0],
[ -4*f1023 - 4*f1122 + 16*f1123 - 4*f1124 - 4*f1223 - 4416.0],
[ -4*f1123 - 4*f1222 + 16*f1223 - 4*f1224 - 4*f1323 + 1928.0],
[  -4*f1223 - 4*f1322 + 16*f1323 - 4*f1324 - 4*f1423 - 712.0],
[ -4*f1323 - 4*f1422 + 16*f1423 - 4*f1424 - 4*f1523 - 4052.0],
[ -4*f1423 - 4*f1522 + 16*f1523 - 4*f1524 - 4*f1623 + 8308.0],
[ -4*f1523 - 4*f1622 + 16*f1623 - 4*f1624 - 4*f1723 - 9532.0],
[ -4*f1623 - 4*f1722 + 16*f1723 - 4*f1724 - 4*f1823 + 2440.0],
[ -4*f1723 - 4*f1822 + 16*f1823 - 4*f1824 - 4*f1923 + 7080.0],
[ -4*f1823 - 4*f1922 + 16*f1923 - 4*f1924 - 4*f2023 - 6792.0],
[ -4*f1923 - 4*f2022 + 16*f2023 - 4*f2024 - 4*f2123 + 1456.0],
[                     -4*f2023 - 4*f2122 + 12*f2123 - 2260.0],
[                       -4*f1024 - 4*f923 + 12*f924 - 4888.0],
[            -4*f1023 + 14*f1024 - 4*f1124 - 4*f924 - 2150.0],
[ -4*f1024 - 4*f1123 + 16*f1124 - 4*f1125 - 4*f1224 + 3468.0],
[  -4*f1124 - 4*f1223 + 16*f1224 - 4*f1225 - 4*f1324 + 104.0],
[  -4*f1224 - 4*f1323 + 16*f1324 - 4*f1325 - 4*f1424 + 832.0],
[ -4*f1324 - 4*f1423 + 16*f1424 - 4*f1425 - 4*f1524 + 9792.0],
[-4*f1424 - 4*f1523 + 16*f1524 - 4*f1525 - 4*f1624 - 10372.0],
[ -4*f1524 - 4*f1623 + 16*f1624 - 4*f1625 - 4*f1724 + 4000.0],
[ -4*f1624 - 4*f1723 + 16*f1724 - 4*f1725 - 4*f1824 + 4072.0],
[ -4*f1724 - 4*f1823 + 16*f1824 - 4*f1825 - 4*f1924 - 7752.0],
[           -4*f1824 - 4*f1923 + 14*f1924 - 4*f2024 + 6116.0],
[                     -4*f1924 - 4*f2023 + 12*f2024 - 4400.0],
[                     -4*f1124 + 12*f1125 - 4*f1225 - 8018.0],
[            -4*f1125 - 4*f1224 + 14*f1225 - 4*f1325 + 616.0],
[           -4*f1225 - 4*f1324 + 14*f1325 - 4*f1425 - 5872.0],
[           -4*f1325 - 4*f1424 + 14*f1425 - 4*f1525 - 3732.0],
[           -4*f1425 - 4*f1524 + 14*f1525 - 4*f1625 + 7208.0],
[           -4*f1525 - 4*f1624 + 14*f1625 - 4*f1725 - 7150.0],
[           -4*f1625 - 4*f1724 + 14*f1725 - 4*f1825 - 2676.0],
[                      -4*f1725 - 4*f1824 + 12*f1825 + 594.0]])

In [32]:

    

sympy.diff(eq3, f3[12][5]).subs(symbolic_replacements3)


    

    
Out[32]:





-2*c124 - 4*f115 + 14*f125 - 4*f126 - 4*f135 + 4*g115 + 2*g124 - 14*g125 + 4*g126 + 4*g135

In [ ]:

    




    

In [49]:

    

display_matrix(mat3.subs(symbolic_replacements3))


    

In [33]:

    

mat3s = mat3.subs(replacements3)


    

In [24]:

    

mat3s


    

    
Out[24]:





Matrix([
[                                          -4*f116 + 3480.0],
[                                          -4*f126 + 2786.0],
[                                          -4*f136 + 3458.0],
[                                          -4*f146 + 3446.0],
[                                          -4*f156 - 1592.0],
[                                          -4*f166 + 2488.0],
[                                          -4*f176 + 3922.0],
[                                           -4*f186 + 418.0],
[                                            -4*f97 + 472.0],
[                                 -4*f107 - 4*f116 + 9340.0],
[                        16*f116 - 4*f117 - 4*f126 - 7582.0],
[              -4*f116 + 16*f126 - 4*f127 - 4*f136 - 4774.0],
[              -4*f126 + 16*f136 - 4*f137 - 4*f146 - 1420.0],
[              -4*f136 + 16*f146 - 4*f147 - 4*f156 - 1582.0],
[              -4*f146 + 16*f156 - 4*f157 - 4*f166 - 2808.0],
[               -4*f156 + 16*f166 - 4*f167 - 4*f176 + 702.0],
[              -4*f166 + 16*f176 - 4*f177 - 4*f186 - 3410.0],
[                       -4*f176 + 16*f186 - 4*f187 - 2554.0],
[                                 -4*f186 - 4*f197 + 2120.0],
[                                           -4*f207 + 398.0],
[                                   -4*f88 - 4*f97 + 3946.0],
[                         -4*f107 + 16*f97 - 4*f98 - 3108.0],
[                16*f107 - 4*f108 - 4*f117 - 4*f97 - 5500.0],
[     -4*f107 - 4*f116 + 16*f117 - 4*f118 - 4*f127 + 2464.0],
[      -4*f117 - 4*f126 + 16*f127 - 4*f128 - 4*f137 + 318.0],
[     -4*f127 - 4*f136 + 16*f137 - 4*f138 - 4*f147 + 1586.0],
[     -4*f137 - 4*f146 + 16*f147 - 4*f148 - 4*f157 - 3554.0],
[     -4*f147 - 4*f156 + 16*f157 - 4*f158 - 4*f167 + 2198.0],
[     -4*f157 - 4*f166 + 16*f167 - 4*f168 - 4*f177 - 1000.0],
[     -4*f167 - 4*f176 + 16*f177 - 4*f178 - 4*f187 - 3466.0],
[     -4*f177 - 4*f186 + 16*f187 - 4*f188 - 4*f197 + 2494.0],
[               -4*f187 + 16*f197 - 4*f198 - 4*f207 - 858.0],
[                       -4*f197 + 16*f207 - 4*f208 - 2474.0],
[                                 -4*f207 - 4*f218 - 3024.0],
[                                   -4*f79 - 4*f88 + 7668.0],
[                           16*f88 - 4*f89 - 4*f98 - 6736.0],
[           -4*f108 - 4*f88 - 4*f97 + 16*f98 - 4*f99 + 40.0],
[       -4*f107 + 16*f108 - 4*f109 - 4*f118 - 4*f98 + 350.0],
[     -4*f108 - 4*f117 + 16*f118 - 4*f119 - 4*f128 - 1478.0],
[     -4*f118 - 4*f127 + 16*f128 - 4*f129 - 4*f138 + 1652.0],
[     -4*f128 - 4*f137 + 16*f138 - 4*f139 - 4*f148 - 4812.0],
[     -4*f138 - 4*f147 + 16*f148 - 4*f149 - 4*f158 + 4366.0],
[     -4*f148 - 4*f157 + 16*f158 - 4*f159 - 4*f168 + 1874.0],
[     -4*f158 - 4*f167 + 16*f168 - 4*f169 - 4*f178 - 2612.0],
[     -4*f168 - 4*f177 + 16*f178 - 4*f179 - 4*f188 + 2558.0],
[     -4*f178 - 4*f187 + 16*f188 - 4*f189 - 4*f198 - 1966.0],
[     -4*f188 - 4*f197 + 16*f198 - 4*f199 - 4*f208 + 1304.0],
[      -4*f198 - 4*f207 + 16*f208 - 4*f209 - 4*f218 - 580.0],
[                        -4*f208 + 16*f218 - 4*f219 + 518.0],
[                                          -4*f218 - 3610.0],
[                                  -4*f610 - 4*f79 + 7250.0],
[                         -4*f710 + 16*f79 - 4*f89 - 9226.0],
[         -4*f79 - 4*f810 - 4*f88 + 16*f89 - 4*f99 + 2282.0],
[         -4*f109 - 4*f89 - 4*f910 - 4*f98 + 16*f99 - 268.0],
[      -4*f1010 - 4*f108 + 16*f109 - 4*f119 - 4*f99 + 322.0],
[    -4*f109 - 4*f1110 - 4*f118 + 16*f119 - 4*f129 + 1634.0],
[    -4*f119 - 4*f1210 - 4*f128 + 16*f129 - 4*f139 - 2622.0],
[    -4*f129 - 4*f1310 - 4*f138 + 16*f139 - 4*f149 + 5136.0],
[    -4*f139 - 4*f1410 - 4*f148 + 16*f149 - 4*f159 - 1132.0],
[    -4*f149 - 4*f1510 - 4*f158 + 16*f159 - 4*f169 - 6516.0],
[    -4*f159 - 4*f1610 - 4*f168 + 16*f169 - 4*f179 + 7218.0],
[    -4*f169 - 4*f1710 - 4*f178 + 16*f179 - 4*f189 - 3284.0],
[     -4*f179 - 4*f1810 - 4*f188 + 16*f189 - 4*f199 - 924.0],
[    -4*f189 - 4*f1910 - 4*f198 + 16*f199 - 4*f209 + 1168.0],
[     -4*f199 - 4*f2010 - 4*f208 + 16*f209 - 4*f219 - 838.0],
[             -4*f209 - 4*f2110 - 4*f218 + 16*f219 - 1326.0],
[                                -4*f219 - 4*f2210 + 2224.0],
[                                           -4*f610 + 976.0],
[                        16*f610 - 4*f611 - 4*f710 - 6466.0],
[       -4*f610 + 16*f710 - 4*f711 - 4*f79 - 4*f810 + 746.0],
[        -4*f710 + 16*f810 - 4*f811 - 4*f89 - 4*f910 - 62.0],
[     -4*f1010 - 4*f810 + 16*f910 - 4*f911 - 4*f99 + 1200.0],
[   16*f1010 - 4*f1011 - 4*f109 - 4*f1110 - 4*f910 + 1092.0],
[  -4*f1010 + 16*f1110 - 4*f1111 - 4*f119 - 4*f1210 - 756.0],
[  -4*f1110 + 16*f1210 - 4*f1211 - 4*f129 - 4*f1310 - 836.0],
[ -4*f1210 + 16*f1310 - 4*f1311 - 4*f139 - 4*f1410 - 1228.0],
[  -4*f1310 + 16*f1410 - 4*f1411 - 4*f149 - 4*f1510 - 770.0],
[ -4*f1410 + 16*f1510 - 4*f1511 - 4*f159 - 4*f1610 + 1996.0],
[  -4*f1510 + 16*f1610 - 4*f1611 - 4*f169 - 4*f1710 + 664.0],
[ -4*f1610 + 16*f1710 - 4*f1711 - 4*f179 - 4*f1810 - 3088.0],
[ -4*f1710 + 16*f1810 - 4*f1811 - 4*f189 - 4*f1910 + 3088.0],
[ -4*f1810 + 16*f1910 - 4*f1911 - 4*f199 - 4*f2010 - 2048.0],
[ -4*f1910 + 16*f2010 - 4*f2011 - 4*f209 - 4*f2110 + 2196.0],
[  -4*f2010 + 16*f2110 - 4*f2111 - 4*f219 - 4*f2210 - 548.0],
[                     -4*f2110 + 16*f2210 - 4*f2211 - 508.0],
[                                         -4*f2210 - 3376.0],
[                                          -4*f611 + 4300.0],
[               -4*f610 + 16*f611 - 4*f612 - 4*f711 + 340.0],
[      -4*f611 - 4*f710 + 16*f711 - 4*f712 - 4*f811 + 174.0],
[      -4*f711 - 4*f810 + 16*f811 - 4*f812 - 4*f911 - 164.0],
[    -4*f1011 - 4*f811 - 4*f910 + 16*f911 - 4*f912 - 3392.0],
[  -4*f1010 + 16*f1011 - 4*f1012 - 4*f1111 - 4*f911 - 582.0],
[ -4*f1011 - 4*f1110 + 16*f1111 - 4*f1112 - 4*f1211 - 476.0],
[-4*f1111 - 4*f1210 + 16*f1211 - 4*f1212 - 4*f1311 + 2024.0],
[-4*f1211 - 4*f1310 + 16*f1311 - 4*f1312 - 4*f1411 - 1526.0],
[-4*f1311 - 4*f1410 + 16*f1411 - 4*f1412 - 4*f1511 + 1664.0],
[ -4*f1411 - 4*f1510 + 16*f1511 - 4*f1512 - 4*f1611 - 308.0],
[-4*f1511 - 4*f1610 + 16*f1611 - 4*f1612 - 4*f1711 - 2384.0],
[-4*f1611 - 4*f1710 + 16*f1711 - 4*f1712 - 4*f1811 + 1544.0],
[ -4*f1711 - 4*f1810 + 16*f1811 - 4*f1812 - 4*f1911 + 298.0],
[ -4*f1811 - 4*f1910 + 16*f1911 - 4*f1912 - 4*f2011 + 254.0],
[-4*f1911 - 4*f2010 + 16*f2011 - 4*f2012 - 4*f2111 - 2544.0],
[-4*f2011 - 4*f2110 + 16*f2111 - 4*f2112 - 4*f2211 + 1346.0],
[          -4*f2111 - 4*f2210 + 16*f2211 - 4*f2212 - 1298.0],
[                               -4*f2211 - 4*f2312 + 1618.0],
[                                  -4*f513 - 4*f612 - 232.0],
[                -4*f611 + 16*f612 - 4*f613 - 4*f712 + 76.0],
[     -4*f612 - 4*f711 + 16*f712 - 4*f713 - 4*f812 - 2858.0],
[     -4*f712 - 4*f811 + 16*f812 - 4*f813 - 4*f912 + 1008.0],
[     -4*f1012 - 4*f812 - 4*f911 + 16*f912 - 4*f913 + 712.0],
[  -4*f1011 + 16*f1012 - 4*f1013 - 4*f1112 - 4*f912 - 306.0],
[-4*f1012 - 4*f1111 + 16*f1112 - 4*f1113 - 4*f1212 + 1556.0],
[ -4*f1112 - 4*f1211 + 16*f1212 - 4*f1213 - 4*f1312 - 958.0],
[-4*f1212 - 4*f1311 + 16*f1312 - 4*f1313 - 4*f1412 + 2608.0],
[-4*f1312 - 4*f1411 + 16*f1412 - 4*f1413 - 4*f1512 - 2134.0],
[ -4*f1412 - 4*f1511 + 16*f1512 - 4*f1513 - 4*f1612 + 916.0],
[ -4*f1512 - 4*f1611 + 16*f1612 - 4*f1613 - 4*f1712 + 472.0],
[ -4*f1612 - 4*f1711 + 16*f1712 - 4*f1713 - 4*f1812 + 664.0],
[ -4*f1712 - 4*f1811 + 16*f1812 - 4*f1813 - 4*f1912 - 726.0],
[-4*f1812 - 4*f1911 + 16*f1912 - 4*f1913 - 4*f2012 + 2816.0],
[ -4*f1912 - 4*f2011 + 16*f2012 - 4*f2013 - 4*f2112 - 444.0],
[-4*f2012 - 4*f2111 + 16*f2112 - 4*f2113 - 4*f2212 - 1946.0],
[-4*f2112 - 4*f2211 + 16*f2212 - 4*f2213 - 4*f2312 + 1148.0],
[                    -4*f2212 + 16*f2312 - 4*f2313 - 2412.0],
[                               -4*f2312 - 4*f2413 - 1646.0],
[                                          -4*f513 + 5282.0],
[                        16*f513 - 4*f514 - 4*f613 - 4482.0],
[      -4*f513 - 4*f612 + 16*f613 - 4*f614 - 4*f713 - 292.0],
[       -4*f613 - 4*f712 + 16*f713 - 4*f714 - 4*f813 - 10.0],
[     -4*f713 - 4*f812 + 16*f813 - 4*f814 - 4*f913 + 1854.0],
[     -4*f1013 - 4*f813 - 4*f912 + 16*f913 - 4*f914 + 790.0],
[ -4*f1012 + 16*f1013 - 4*f1014 - 4*f1113 - 4*f913 - 1424.0],
[ -4*f1013 - 4*f1112 + 16*f1113 - 4*f1114 - 4*f1213 + 426.0],
[-4*f1113 - 4*f1212 + 16*f1213 - 4*f1214 - 4*f1313 - 1008.0],
[ -4*f1213 - 4*f1312 + 16*f1313 - 4*f1314 - 4*f1413 + 156.0],
[ -4*f1313 - 4*f1412 + 16*f1413 - 4*f1414 - 4*f1513 - 432.0],
[ -4*f1413 - 4*f1512 + 16*f1513 - 4*f1514 - 4*f1613 + 190.0],
[ -4*f1513 - 4*f1612 + 16*f1613 - 4*f1614 - 4*f1713 - 274.0],
[ -4*f1613 - 4*f1712 + 16*f1713 - 4*f1714 - 4*f1813 + 414.0],
[-4*f1713 - 4*f1812 + 16*f1813 - 4*f1814 - 4*f1913 - 1000.0],
[-4*f1813 - 4*f1912 + 16*f1913 - 4*f1914 - 4*f2013 - 3370.0],
[-4*f1913 - 4*f2012 + 16*f2013 - 4*f2014 - 4*f2113 + 2580.0],
[-4*f2013 - 4*f2112 + 16*f2113 - 4*f2114 - 4*f2213 + 1716.0],
[-4*f2113 - 4*f2212 + 16*f2213 - 4*f2214 - 4*f2313 - 1052.0],
[-4*f2213 - 4*f2312 + 16*f2313 - 4*f2314 - 4*f2413 + 1170.0],
[                    -4*f2313 + 16*f2413 - 4*f2414 - 1574.0],
[                                         -4*f2413 - 2650.0],
[                                          -4*f514 - 2940.0],
[              -4*f513 + 16*f514 - 4*f515 - 4*f614 + 1550.0],
[     -4*f514 - 4*f613 + 16*f614 - 4*f615 - 4*f714 - 2876.0],
[     -4*f614 - 4*f713 + 16*f714 - 4*f715 - 4*f814 + 1812.0],
[     -4*f714 - 4*f813 + 16*f814 - 4*f815 - 4*f914 - 2200.0],
[     -4*f1014 - 4*f814 - 4*f913 + 16*f914 - 4*f915 - 824.0],
[ -4*f1013 + 16*f1014 - 4*f1015 - 4*f1114 - 4*f914 + 2402.0],
[ -4*f1014 - 4*f1113 + 16*f1114 - 4*f1115 - 4*f1214 - 912.0],
[ -4*f1114 - 4*f1213 + 16*f1214 - 4*f1215 - 4*f1314 + 166.0],
[ -4*f1214 - 4*f1313 + 16*f1314 - 4*f1315 - 4*f1414 - 416.0],
[-4*f1314 - 4*f1413 + 16*f1414 - 4*f1415 - 4*f1514 + 1938.0],
[-4*f1414 - 4*f1513 + 16*f1514 - 4*f1515 - 4*f1614 - 1544.0],
[-4*f1514 - 4*f1613 + 16*f1614 - 4*f1615 - 4*f1714 - 3372.0],
[-4*f1614 - 4*f1713 + 16*f1714 - 4*f1715 - 4*f1814 + 2990.0],
[-4*f1714 - 4*f1813 + 16*f1814 - 4*f1815 - 4*f1914 - 1380.0],
[-4*f1814 - 4*f1913 + 16*f1914 - 4*f1915 - 4*f2014 + 1770.0],
[ -4*f1914 - 4*f2013 + 16*f2014 - 4*f2015 - 4*f2114 - 562.0],
[ -4*f2014 - 4*f2113 + 16*f2114 - 4*f2115 - 4*f2214 - 416.0],
[  -4*f2114 - 4*f2213 + 16*f2214 - 4*f2215 - 4*f2314 - 76.0],
[ -4*f2214 - 4*f2313 + 16*f2314 - 4*f2315 - 4*f2414 - 384.0],
[           -4*f2314 - 4*f2413 + 16*f2414 - 4*f2415 - 522.0],
[                                          -4*f2414 + 182.0],
[                                          -4*f515 + 4118.0],
[              -4*f514 + 16*f515 - 4*f516 - 4*f615 - 4438.0],
[     -4*f515 - 4*f614 + 16*f615 - 4*f616 - 4*f715 + 3176.0],
[     -4*f615 - 4*f714 + 16*f715 - 4*f716 - 4*f815 - 2618.0],
[     -4*f715 - 4*f814 + 16*f815 - 4*f816 - 4*f915 + 2182.0],
[     -4*f1015 - 4*f815 - 4*f914 + 16*f915 - 4*f916 - 766.0],
[ -4*f1014 + 16*f1015 - 4*f1016 - 4*f1115 - 4*f915 - 1654.0],
[-4*f1015 - 4*f1114 + 16*f1115 - 4*f1116 - 4*f1215 + 1098.0],
[  -4*f1115 - 4*f1214 + 16*f1215 - 4*f1216 - 4*f1315 + 66.0],
[-4*f1215 - 4*f1314 + 16*f1315 - 4*f1316 - 4*f1415 - 2274.0],
[ -4*f1315 - 4*f1414 + 16*f1415 - 4*f1416 - 4*f1515 + 498.0],
[ -4*f1415 - 4*f1514 + 16*f1515 - 4*f1516 - 4*f1615 + 820.0],
[-4*f1515 - 4*f1614 + 16*f1615 - 4*f1616 - 4*f1715 + 4056.0],
[-4*f1615 - 4*f1714 + 16*f1715 - 4*f1716 - 4*f1815 + 1790.0],
[-4*f1715 - 4*f1814 + 16*f1815 - 4*f1816 - 4*f1915 - 6860.0],
[-4*f1815 - 4*f1914 + 16*f1915 - 4*f1916 - 4*f2015 + 5900.0],
[-4*f1915 - 4*f2014 + 16*f2015 - 4*f2016 - 4*f2115 - 5224.0],
[ -4*f2015 - 4*f2114 + 16*f2115 - 4*f2116 - 4*f2215 + 980.0],
[ -4*f2115 - 4*f2214 + 16*f2215 - 4*f2216 - 4*f2315 + 498.0],
[  -4*f2215 - 4*f2314 + 16*f2315 - 4*f2316 - 4*f2415 + 68.0],
[          -4*f2315 - 4*f2414 + 16*f2415 - 4*f2416 - 1042.0],
[                                         -4*f2415 + 2170.0],
[                                          -4*f516 - 3086.0],
[                        -4*f515 + 16*f516 - 4*f616 + 986.0],
[     -4*f516 - 4*f615 + 16*f616 - 4*f617 - 4*f716 - 3234.0],
[      -4*f616 - 4*f715 + 16*f716 - 4*f717 - 4*f816 + 710.0],
[      -4*f716 - 4*f815 + 16*f816 - 4*f817 - 4*f916 - 660.0],
[     -4*f1016 - 4*f816 - 4*f915 + 16*f916 - 4*f917 + 470.0],
[ -4*f1015 + 16*f1016 - 4*f1017 - 4*f1116 - 4*f916 + 1428.0],
[ -4*f1016 - 4*f1115 + 16*f1116 - 4*f1117 - 4*f1216 - 744.0],
[ -4*f1116 - 4*f1215 + 16*f1216 - 4*f1217 - 4*f1316 - 256.0],
[-4*f1216 - 4*f1315 + 16*f1316 - 4*f1317 - 4*f1416 + 3408.0],
[-4*f1316 - 4*f1415 + 16*f1416 - 4*f1417 - 4*f1516 + 1318.0],
[-4*f1416 - 4*f1515 + 16*f1516 - 4*f1517 - 4*f1616 - 5076.0],
[  -4*f1516 - 4*f1615 + 16*f1616 - 4*f1617 - 4*f1716 - 44.0],
[ -4*f1616 - 4*f1715 + 16*f1716 - 4*f1717 - 4*f1816 - 372.0],
[-4*f1716 - 4*f1815 + 16*f1816 - 4*f1817 - 4*f1916 + 3004.0],
[-4*f1816 - 4*f1915 + 16*f1916 - 4*f1917 - 4*f2016 - 2852.0],
[-4*f1916 - 4*f2015 + 16*f2016 - 4*f2017 - 4*f2116 + 3126.0],
[ -4*f2016 - 4*f2115 + 16*f2116 - 4*f2117 - 4*f2216 + 128.0],
[ -4*f2116 - 4*f2215 + 16*f2216 - 4*f2217 - 4*f2316 - 506.0],
[ -4*f2216 - 4*f2315 + 16*f2316 - 4*f2317 - 4*f2416 - 104.0],
[          -4*f2316 - 4*f2415 + 16*f2416 - 4*f2417 - 2118.0],
[                                         -4*f2416 + 1206.0],
[                                 -4*f516 - 4*f617 - 3098.0],
[                        -4*f616 + 16*f617 - 4*f717 + 750.0],
[      -4*f617 - 4*f716 + 16*f717 - 4*f718 - 4*f817 + 580.0],
[      -4*f717 - 4*f816 + 16*f817 - 4*f818 - 4*f917 - 346.0],
[     -4*f1017 - 4*f817 - 4*f916 + 16*f917 - 4*f918 - 964.0],
[ -4*f1016 + 16*f1017 - 4*f1018 - 4*f1117 - 4*f917 + 1450.0],
[ -4*f1017 - 4*f1116 + 16*f1117 - 4*f1118 - 4*f1217 - 354.0],
[-4*f1117 - 4*f1216 + 16*f1217 - 4*f1218 - 4*f1317 - 1532.0],
[-4*f1217 - 4*f1316 + 16*f1317 - 4*f1318 - 4*f1417 - 2258.0],
[-4*f1317 - 4*f1416 + 16*f1417 - 4*f1418 - 4*f1517 - 1458.0],
[-4*f1417 - 4*f1516 + 16*f1517 - 4*f1518 - 4*f1617 + 4624.0],
[ -4*f1517 - 4*f1616 + 16*f1617 - 4*f1618 - 4*f1717 + 218.0],
[-4*f1617 - 4*f1716 + 16*f1717 - 4*f1718 - 4*f1817 - 2402.0],
[-4*f1717 - 4*f1816 + 16*f1817 - 4*f1818 - 4*f1917 + 1410.0],
[-4*f1817 - 4*f1916 + 16*f1917 - 4*f1918 - 4*f2017 - 2002.0],
[ -4*f1917 - 4*f2016 + 16*f2017 - 4*f2018 - 4*f2117 - 338.0],
[-4*f2017 - 4*f2116 + 16*f2117 - 4*f2118 - 4*f2217 + 1452.0],
[ -4*f2117 - 4*f2216 + 16*f2217 - 4*f2218 - 4*f2317 - 556.0],
[ -4*f2217 - 4*f2316 + 16*f2317 - 4*f2318 - 4*f2417 - 172.0],
[                    -4*f2317 - 4*f2416 + 16*f2417 - 3632.0],
[                                         -4*f2417 + 2352.0],
[                                 -4*f617 - 4*f718 + 1904.0],
[              -4*f717 + 16*f718 - 4*f719 - 4*f818 - 2182.0],
[     -4*f718 - 4*f817 + 16*f818 - 4*f819 - 4*f918 - 1266.0],
[    -4*f1018 - 4*f818 - 4*f917 + 16*f918 - 4*f919 + 2026.0],
[ -4*f1017 + 16*f1018 - 4*f1019 - 4*f1118 - 4*f918 - 2752.0],
[-4*f1018 - 4*f1117 + 16*f1118 - 4*f1119 - 4*f1218 + 2024.0],
[-4*f1118 - 4*f1217 + 16*f1218 - 4*f1219 - 4*f1318 + 1328.0],
[ -4*f1218 - 4*f1317 + 16*f1318 - 4*f1319 - 4*f1418 - 288.0],
[ -4*f1318 - 4*f1417 + 16*f1418 - 4*f1419 - 4*f1518 + 648.0],
[-4*f1418 - 4*f1517 + 16*f1518 - 4*f1519 - 4*f1618 - 1346.0],
[-4*f1518 - 4*f1617 + 16*f1618 - 4*f1619 - 4*f1718 - 2636.0],
[-4*f1618 - 4*f1717 + 16*f1718 - 4*f1719 - 4*f1818 + 3088.0],
[-4*f1718 - 4*f1817 + 16*f1818 - 4*f1819 - 4*f1918 + 1050.0],
[ -4*f1818 - 4*f1917 + 16*f1918 - 4*f1919 - 4*f2018 - 444.0],
[-4*f1918 - 4*f2017 + 16*f2018 - 4*f2019 - 4*f2118 + 1384.0],
[-4*f2018 - 4*f2117 + 16*f2118 - 4*f2119 - 4*f2218 - 2548.0],
[ -4*f2118 - 4*f2217 + 16*f2218 - 4*f2219 - 4*f2318 + 338.0],
[           -4*f2218 - 4*f2317 + 16*f2318 - 4*f2319 - 744.0],
[                               -4*f2318 - 4*f2417 + 3594.0],
[                                           -4*f719 - 386.0],
[                       -4*f718 + 16*f719 - 4*f819 - 2032.0],
[     -4*f719 - 4*f818 + 16*f819 - 4*f820 - 4*f919 - 2104.0],
[     -4*f1019 - 4*f819 - 4*f918 + 16*f919 - 4*f920 + 614.0],
[  -4*f1018 + 16*f1019 - 4*f1020 - 4*f1119 - 4*f919 - 472.0],
[ -4*f1019 - 4*f1118 + 16*f1119 - 4*f1120 - 4*f1219 - 740.0],
[-4*f1119 - 4*f1218 + 16*f1219 - 4*f1220 - 4*f1319 - 2714.0],
[-4*f1219 - 4*f1318 + 16*f1319 - 4*f1320 - 4*f1419 + 5066.0],
[-4*f1319 - 4*f1418 + 16*f1419 - 4*f1420 - 4*f1519 - 2706.0],
[ -4*f1419 - 4*f1518 + 16*f1519 - 4*f1520 - 4*f1619 + 386.0],
[-4*f1519 - 4*f1618 + 16*f1619 - 4*f1620 - 4*f1719 + 4204.0],
[-4*f1619 - 4*f1718 + 16*f1719 - 4*f1720 - 4*f1819 - 3986.0],
[-4*f1719 - 4*f1818 + 16*f1819 - 4*f1820 - 4*f1919 - 1842.0],
[-4*f1819 - 4*f1918 + 16*f1919 - 4*f1920 - 4*f2019 + 3610.0],
[-4*f1919 - 4*f2018 + 16*f2019 - 4*f2020 - 4*f2119 - 1564.0],
[-4*f2019 - 4*f2118 + 16*f2119 - 4*f2120 - 4*f2219 - 1042.0],
[-4*f2119 - 4*f2218 + 16*f2219 - 4*f2220 - 4*f2319 + 2692.0],
[          -4*f2219 - 4*f2318 + 16*f2319 - 4*f2320 - 3634.0],
[                                         -4*f2319 + 2252.0],
[                                 -4*f719 - 4*f820 + 5842.0],
[               -4*f819 + 16*f820 - 4*f821 - 4*f920 + 798.0],
[    -4*f1020 - 4*f820 - 4*f919 + 16*f920 - 4*f921 - 2438.0],
[ -4*f1019 + 16*f1020 - 4*f1021 - 4*f1120 - 4*f920 + 3146.0],
[ -4*f1020 - 4*f1119 + 16*f1120 - 4*f1121 - 4*f1220 - 650.0],
[-4*f1120 - 4*f1219 + 16*f1220 - 4*f1221 - 4*f1320 + 3240.0],
[-4*f1220 - 4*f1319 + 16*f1320 - 4*f1321 - 4*f1420 - 7902.0],
[-4*f1320 - 4*f1419 + 16*f1420 - 4*f1421 - 4*f1520 + 8032.0],
[-4*f1420 - 4*f1519 + 16*f1520 - 4*f1521 - 4*f1620 - 4796.0],
[-4*f1520 - 4*f1619 + 16*f1620 - 4*f1621 - 4*f1720 - 2354.0],
[-4*f1620 - 4*f1719 + 16*f1720 - 4*f1721 - 4*f1820 + 2390.0],
[-4*f1720 - 4*f1819 + 16*f1820 - 4*f1821 - 4*f1920 + 1330.0],
[-4*f1820 - 4*f1919 + 16*f1920 - 4*f1921 - 4*f2020 - 4554.0],
[-4*f1920 - 4*f2019 + 16*f2020 - 4*f2021 - 4*f2120 + 2516.0],
[ -4*f2020 - 4*f2119 + 16*f2120 - 4*f2121 - 4*f2220 + 448.0],
[-4*f2120 - 4*f2219 + 16*f2220 - 4*f2221 - 4*f2320 - 1294.0],
[                    -4*f2220 - 4*f2319 + 16*f2320 - 3108.0],
[                                         -4*f2320 + 1302.0],
[                                          -4*f821 - 8212.0],
[                       -4*f820 + 16*f821 - 4*f921 + 1646.0],
[    -4*f1021 - 4*f821 - 4*f920 + 16*f921 - 4*f922 - 3390.0],
[   -4*f1020 + 16*f1021 - 4*f1022 - 4*f1121 - 4*f921 + 78.0],
[-4*f1021 - 4*f1120 + 16*f1121 - 4*f1122 - 4*f1221 - 1880.0],
[-4*f1121 - 4*f1220 + 16*f1221 - 4*f1222 - 4*f1321 + 1560.0],
[-4*f1221 - 4*f1320 + 16*f1321 - 4*f1322 - 4*f1421 + 1068.0],
[-4*f1321 - 4*f1420 + 16*f1421 - 4*f1422 - 4*f1521 - 3592.0],
[-4*f1421 - 4*f1520 + 16*f1521 - 4*f1522 - 4*f1621 + 6112.0],
[-4*f1521 - 4*f1620 + 16*f1621 - 4*f1622 - 4*f1721 - 1286.0],
[ -4*f1621 - 4*f1720 + 16*f1721 - 4*f1722 - 4*f1821 + 182.0],
[-4*f1721 - 4*f1820 + 16*f1821 - 4*f1822 - 4*f1921 - 1050.0],
[-4*f1821 - 4*f1920 + 16*f1921 - 4*f1922 - 4*f2021 + 3666.0],
[-4*f1921 - 4*f2020 + 16*f2021 - 4*f2022 - 4*f2121 - 2086.0],
[ -4*f2021 - 4*f2120 + 16*f2121 - 4*f2122 - 4*f2221 + 752.0],
[                    -4*f2121 - 4*f2220 + 16*f2221 - 3404.0],
[                               -4*f2221 - 4*f2320 + 3666.0],
[                                 -4*f821 - 4*f922 - 1016.0],
[                      -4*f1022 - 4*f921 + 16*f922 + 2076.0],
[ -4*f1021 + 16*f1022 - 4*f1023 - 4*f1122 - 4*f922 - 1452.0],
[-4*f1022 - 4*f1121 + 16*f1122 - 4*f1123 - 4*f1222 + 1320.0],
[-4*f1122 - 4*f1221 + 16*f1222 - 4*f1223 - 4*f1322 - 2492.0],
[-4*f1222 - 4*f1321 + 16*f1322 - 4*f1323 - 4*f1422 + 1572.0],
[ -4*f1322 - 4*f1421 + 16*f1422 - 4*f1423 - 4*f1522 - 646.0],
[-4*f1422 - 4*f1521 + 16*f1522 - 4*f1523 - 4*f1622 - 2438.0],
[-4*f1522 - 4*f1621 + 16*f1622 - 4*f1623 - 4*f1722 + 1740.0],
[ -4*f1622 - 4*f1721 + 16*f1722 - 4*f1723 - 4*f1822 + 494.0],
[-4*f1722 - 4*f1821 + 16*f1822 - 4*f1823 - 4*f1922 - 2904.0],
[ -4*f1822 - 4*f1921 + 16*f1922 - 4*f1923 - 4*f2022 + 868.0],
[ -4*f1922 - 4*f2021 + 16*f2022 - 4*f2023 - 4*f2122 - 344.0],
[                    -4*f2022 - 4*f2121 + 16*f2122 - 2524.0],
[                               -4*f2122 - 4*f2221 + 2280.0],
[                                -4*f1023 - 4*f922 - 1514.0],
[                      -4*f1022 + 16*f1023 - 4*f1123 + 68.0],
[-4*f1023 - 4*f1122 + 16*f1123 - 4*f1124 - 4*f1223 - 2208.0],
[ -4*f1123 - 4*f1222 + 16*f1223 - 4*f1224 - 4*f1323 + 964.0],
[ -4*f1223 - 4*f1322 + 16*f1323 - 4*f1324 - 4*f1423 - 356.0],
[-4*f1323 - 4*f1422 + 16*f1423 - 4*f1424 - 4*f1523 - 2026.0],
[-4*f1423 - 4*f1522 + 16*f1523 - 4*f1524 - 4*f1623 + 4154.0],
[-4*f1523 - 4*f1622 + 16*f1623 - 4*f1624 - 4*f1723 - 4766.0],
[-4*f1623 - 4*f1722 + 16*f1723 - 4*f1724 - 4*f1823 + 1220.0],
[-4*f1723 - 4*f1822 + 16*f1823 - 4*f1824 - 4*f1923 + 3540.0],
[          -4*f1823 - 4*f1922 + 16*f1923 - 4*f2023 - 5296.0],
[                    -4*f1923 - 4*f2022 + 16*f2023 - 3124.0],
[                               -4*f2023 - 4*f2122 + 3822.0],
[                                                   -1834.0],
[                               -4*f1023 - 4*f1124 + 4086.0],
[                    -4*f1123 + 16*f1124 - 4*f1224 - 2958.0],
[          -4*f1124 - 4*f1223 + 16*f1224 - 4*f1324 - 2148.0],
[          -4*f1224 - 4*f1323 + 16*f1324 - 4*f1424 - 2028.0],
[          -4*f1324 - 4*f1423 + 16*f1424 - 4*f1524 + 2820.0],
[          -4*f1424 - 4*f1523 + 16*f1524 - 4*f1624 - 6430.0],
[           -4*f1524 - 4*f1623 + 16*f1624 - 4*f1724 - 392.0],
[            -4*f1624 - 4*f1723 + 16*f1724 - 4*f1824 - 56.0],
[                    -4*f1724 - 4*f1823 + 16*f1824 - 7972.0],
[                               -4*f1824 - 4*f1923 + 7812.0],
[                                          -4*f2023 - 198.0],
[                                          -4*f1124 - 758.0],
[                                         -4*f1224 + 2814.0],
[                                         -4*f1324 + 1052.0],
[                                          -4*f1424 + 698.0],
[                                         -4*f1524 + 3444.0],
[                                          -4*f1624 + 698.0],
[                                         -4*f1724 + 1026.0],
[                                         -4*f1824 + 3304.0]])

In [26]:

    

vec3


    

    
Out[26]:





[f115,
 f125,
 f135,
 f145,
 f155,
 f165,
 f175,
 f185,
 f96,
 f106,
 f116,
 f126,
 f136,
 f146,
 f156,
 f166,
 f176,
 f186,
 f196,
 f206,
 f87,
 f97,
 f107,
 f117,
 f127,
 f137,
 f147,
 f157,
 f167,
 f177,
 f187,
 f197,
 f207,
 f217,
 f78,
 f88,
 f98,
 f108,
 f118,
 f128,
 f138,
 f148,
 f158,
 f168,
 f178,
 f188,
 f198,
 f208,
 f218,
 f228,
 f69,
 f79,
 f89,
 f99,
 f109,
 f119,
 f129,
 f139,
 f149,
 f159,
 f169,
 f179,
 f189,
 f199,
 f209,
 f219,
 f229,
 f510,
 f610,
 f710,
 f810,
 f910,
 f1010,
 f1110,
 f1210,
 f1310,
 f1410,
 f1510,
 f1610,
 f1710,
 f1810,
 f1910,
 f2010,
 f2110,
 f2210,
 f2310,
 f511,
 f611,
 f711,
 f811,
 f911,
 f1011,
 f1111,
 f1211,
 f1311,
 f1411,
 f1511,
 f1611,
 f1711,
 f1811,
 f1911,
 f2011,
 f2111,
 f2211,
 f2311,
 f512,
 f612,
 f712,
 f812,
 f912,
 f1012,
 f1112,
 f1212,
 f1312,
 f1412,
 f1512,
 f1612,
 f1712,
 f1812,
 f1912,
 f2012,
 f2112,
 f2212,
 f2312,
 f2412,
 f413,
 f513,
 f613,
 f713,
 f813,
 f913,
 f1013,
 f1113,
 f1213,
 f1313,
 f1413,
 f1513,
 f1613,
 f1713,
 f1813,
 f1913,
 f2013,
 f2113,
 f2213,
 f2313,
 f2413,
 f2513,
 f414,
 f514,
 f614,
 f714,
 f814,
 f914,
 f1014,
 f1114,
 f1214,
 f1314,
 f1414,
 f1514,
 f1614,
 f1714,
 f1814,
 f1914,
 f2014,
 f2114,
 f2214,
 f2314,
 f2414,
 f2514,
 f415,
 f515,
 f615,
 f715,
 f815,
 f915,
 f1015,
 f1115,
 f1215,
 f1315,
 f1415,
 f1515,
 f1615,
 f1715,
 f1815,
 f1915,
 f2015,
 f2115,
 f2215,
 f2315,
 f2415,
 f2515,
 f416,
 f516,
 f616,
 f716,
 f816,
 f916,
 f1016,
 f1116,
 f1216,
 f1316,
 f1416,
 f1516,
 f1616,
 f1716,
 f1816,
 f1916,
 f2016,
 f2116,
 f2216,
 f2316,
 f2416,
 f2516,
 f517,
 f617,
 f717,
 f817,
 f917,
 f1017,
 f1117,
 f1217,
 f1317,
 f1417,
 f1517,
 f1617,
 f1717,
 f1817,
 f1917,
 f2017,
 f2117,
 f2217,
 f2317,
 f2417,
 f2517,
 f618,
 f718,
 f818,
 f918,
 f1018,
 f1118,
 f1218,
 f1318,
 f1418,
 f1518,
 f1618,
 f1718,
 f1818,
 f1918,
 f2018,
 f2118,
 f2218,
 f2318,
 f2418,
 f619,
 f719,
 f819,
 f919,
 f1019,
 f1119,
 f1219,
 f1319,
 f1419,
 f1519,
 f1619,
 f1719,
 f1819,
 f1919,
 f2019,
 f2119,
 f2219,
 f2319,
 f2419,
 f720,
 f820,
 f920,
 f1020,
 f1120,
 f1220,
 f1320,
 f1420,
 f1520,
 f1620,
 f1720,
 f1820,
 f1920,
 f2020,
 f2120,
 f2220,
 f2320,
 f2420,
 f721,
 f821,
 f921,
 f1021,
 f1121,
 f1221,
 f1321,
 f1421,
 f1521,
 f1621,
 f1721,
 f1821,
 f1921,
 f2021,
 f2121,
 f2221,
 f2321,
 f822,
 f922,
 f1022,
 f1122,
 f1222,
 f1322,
 f1422,
 f1522,
 f1622,
 f1722,
 f1822,
 f1922,
 f2022,
 f2122,
 f2222,
 f923,
 f1023,
 f1123,
 f1223,
 f1323,
 f1423,
 f1523,
 f1623,
 f1723,
 f1823,
 f1923,
 f2023,
 f2123,
 f924,
 f1024,
 f1124,
 f1224,
 f1324,
 f1424,
 f1524,
 f1624,
 f1724,
 f1824,
 f1924,
 f2024,
 f1125,
 f1225,
 f1325,
 f1425,
 f1525,
 f1625,
 f1725,
 f1825]

In [29]:

    

# x1 = sympy.solve(mat1.subs(replacements1))
# x2 = sympy.solve(mat2.subs(replacements2))
x3 = sympy.solve(mat3s)


    

In [30]:

    

x3


    

    
Out[30]:





[]

In [21]:

    

def apply_blended_pixels(x, reference, f):
    img = np.copy(reference)
    
    height, width = reference.shape
    
    for j in range(height):
        for i in range(width):
            if x.has_key(f[i][j]):
                img[j][i] = x[f[i][j]]
    
    return img


    

In [22]:

    

x3


    

    
Out[22]:





[]

In [52]:

    

# img1 = apply_blended_pixels(x1, reference, f1)
# img2 = apply_blended_pixels(x2, reference, f2)
img3 = apply_blended_pixels(x3, reference, f3)


    

In [ ]:

    

f, (ax1, ax2) = plt.subplots(1, 2, sharex=True, sharey=True, figsize=(12, 4))
ax1.imshow(img1, vmin=300, vmax=750)
ax2.imshow(img2, vmin=300, vmax=750)


    

In [ ]:

    

print img1


    

In [ ]:

    

print img2