Emulating nested for loop with scan is slow

I am trying to convolve over an image with a special (bilateral) filter all while using Tensorflow control ops so I can use automatic differentiation. To accomplish this, I nested scan() functions, but this leaves me with remarkably poor performance---a small image takes more than 5 minutes. Is there a way to speed this operation, and how badly am I using Tensorflow control flow operations?

sess = tf.InteractiveSession()
with sess.as_default():
    def bilateralFilter(image, sigma_space, sigma_range):
        win_size = max(5., 2 * np.ceil(3 * sigma_space) + 1)

        win_ext = int((win_size - 1) / 2)
        height = tf.shape(image)[0].eval()
        width = tf.shape(image)[1].eval()


        spatial_gaussian = []
        for i in range(-win_ext, win_ext + 1):
            for j in range(-win_ext, win_ext + 1):
                spatial_gaussian.append(np.exp(-0.5 *
                (i ** 2 + j ** 2) / sigma_space ** 2))

        # we use "symmetric" as it best approximates "edge" padding
        padded = tf.pad(image, [[win_ext, win_ext],
                 [win_ext, win_ext]], mode='SYMMETRIC')
        global out_image
        out_image = tf.zeros(tf.shape(image))
        global weight
        weight = tf.zeros(tf.shape(image))

        global spatial_index
        spatial_index = 0
        global row
        row = -win_ext
        global col
        col = -win_ext

        elems = list(range(-win_ext, win_ext + 1))

        def cols():
            global out_image
            global weight
            global spatial_index
            global row
            global col
            sub_image = tf.slice(padded,
                      [win_ext + row, win_ext + col], [height, width])
            value = tf.exp(-0.5 *(((image - sub_image) / sigma_range) ** 2))
            * spatial_gaussian[spatial_index]

            out_image += value * sub_image
            weight += value
            spatial_index += 1
            col += 1
            if col % (2*win_ext + 1):
                row += 1    # move to next iteration of larger loop
                col = -win_ext
            retu 0    # we have to retu something

        def rows():
            tf.scan(lambda a, x: cols(), elems)
            retu 0


        tf.scan(lambda a, x: rows(), elems)


        out_image /= weight

        retu out_image

    cat = plt.imread("cat.png")  # grayscale
    cat = tf.reshape(tf.constant(cat), [276, 276])
    cat_blurred = bilateralFilter(cat, 2., 0.25)
    cat_blurred = cat_blurred.eval()
    plt.figure()
    plt.gray()
    plt.imshow(cat_blurred)
    plt.show()

For reference, the above code stems from the following faster code:

def bilateralFilter(image, sigma_space=1, sigma_range=None, win_size=None):

    if sigma_range is None:
        sigma_range = sigma_space
    if win_size is None: win_size = max(5, 2 * int(np.ceil(3*sigma_space)) + 1)

    win_ext = (win_size - 1) / 2
    height = image.shape[0]
    width = image.shape[1]

    # pre-calculate spatial_gaussian
    spatial_gaussian = []
    for i in range(-win_ext, win_ext+1):
        for j in range(-win_ext, win_ext+1):
            spatial_gaussian.append(np.exp(-0.5*(i**2+j**2)/sigma_space**2))

    padded = np.pad(image, win_ext, mode="edge")

    out_image = np.zeros(image.shape)
    weight = np.zeros(image.shape)

    idx = 0
    for row in xrange(-win_ext, 1+win_ext):
        for col in xrange(-win_ext, 1+win_ext):
            slice = padded[win_ext+row:height+win_ext+row,
                                          win_ext+col:width+win_ext+col]
            value = np.exp(-0.5*((image - slice)/sigma_range)**2) 
                    * spatial_gaussian[idx]
            out_image += value*slice
            weight += value
            idx += 1

    out_image /= weight

    retu out_image