CODEWITHSUNDEEP
pythonarraysnumpypixelinterpolation
David Colpey
David Colpey

Reputation: 81

Array/Image Interpolation in python

I have 365 2d numpy arrays for the every day of the year, displaying an image like this:

http://i50.tinypic.com/34i62gw.jpg

I have them all stacked in a 3d numpy array. Pixels with a value that represents cloud i want to get rid of, i want to search through the previous 7 days, or next 7 days (previous 7 layers, next 7 layers) to find a value other than cloud, and then replace the cloud value with other possible values for that pixel (values experienced in the other days/layers for the corresponding pixel).

I am new to python, and am a bit lost.

Any ideas?

Thanks

Upvotes: 3

Views: 917

Answers (3)

David Colpey
David Colpey

Reputation: 81

I solved it by this: 

interpdata = []
j = 0
for i in stack:
    try:
        temp = np.where( stack[j] == 50, stack[j-1], modis[j] )
        temp = np.where( temp == 50, stack[j+1], temp )
        temp = np.where( temp == 50, stack[j-2], temp )
        temp = np.where( temp == 50, stack[j+2], temp )
        temp = np.where( temp == 50, stack[j-3], temp )
        temp = np.where( temp == 50, stack[j+3], temp ) 
        temp = np.where( temp == 50, stack[j-4], temp )
        temp = np.where( temp == 50, stack[j+4], temp )
    except IndexError:
        print 'IndexError Passed'       
        pass
    else:
        pass
    interpdata [j, :, :] = temp
    j = j + 1

Upvotes: 1

Dunes
Dunes

Reputation: 40703

You are essentially trying to write a filter for your array.

First you need to write a function that when given an array of values, with the middle one being the element currently examined, will return some computation of those values. In your case the function will expect to take 1-d array and returns the element nearest to the middle index that is not cloud:

import numpy as np
from scipy.ndimage.filters import generic_filter

_cloud = -1

def findNearestNonCloud(elements):
    middleIndex = len(elements) / 2
    if elements[middleIndex] != _cloud:
        return elements[middleIndex] # middle value is not cloud

    nonCloudIndices, = np.where(elements != _cloud)
    if len(nonCloudIndices) == 0:
        return elements[middleIndex] # all values were cloud

    prevNonCloudIndex = np.where(nonCloudIndices < middleIndex, 
            nonCloudIndices, -1).max()
    nextNonCloudIndex = -np.where(nonCloudIndices > middleIndex, 
            -nonCloudIndices, 1).min()
    # -1 means no non-cloud index

    # pick index closest to middle index    
    if (abs(prevNonCloudIndex - middleIndex) 
            <= abs(nextNonCloudIndex - middleIndex)):
        return elements[prevNonCloudIndex]
    else:
        return elements[nextNonCloudIndex]

Now you need to apply this function to the elements you're interesting. To do this you'll need a mask that denotes which other elements you'll be interested in with regard to a specific element. 

from scipy.ndimage.filters import generic_filter

# creates 5 days worth of a 3x3 plot of land
input = np.ones((5, 3, 3)) * _cloud
input[0,:,:] = 10 # set first "image" to all be 10s
input[4,0,0] = 12 # uppper left corner of fourth image is set to 12
print "input data\n", input, "\n"

mask = (5, 1, 1)
# mask represents looking at the present day, 2 days in the future and 2 days in 
# the past for 5 days in total.

print "result\n", generic_filter(input, findNearestNonCloud, size=mask)
# second and third images should mirror first image, 
# except upper left corner of third image should be 12

Upvotes: 3

mgilson
mgilson

Reputation: 309881

I would think that you could do something like:

data = somehow_get_your_3d_data() #indexed as [day_of_year,y,x]
for i,dat in enumerate(data):
    weeks2 = data[max(i-7,i):min(i+7,len(data)), ... ]
    new_value = get_new_value(weeks2) #get value from weeks2 here somehow
    dat[dat == cloud_value] = new_value

Upvotes: 0

Related Questions