用于3个矩阵的循环

问题描述 投票:0回答:1

我有三个大小相同的2D数组:

import numpy as np
rArray = np.ones((3,3))
gArray = np.ones((3,3))*2
bArray = np.ones((3,3))*3

我想逐个元素地对其进行循环,以便进行一些计算,然后将结果存储在新的数组中:

data = np.zeros(rArray.shape)
for r,g,b in .... :
  res=2*r + g - b
  data[i,j] = res

其中r属于rArray,依此类推。

EDIT1:

# Python program to get average of a list 
def Average(lst): 
    return sum(lst) / len(lst)

# TURBO
from sympy import symbols, solve
import numpy as np

minValue = -140; maxValue = -80

a = 1 / (maxValue - minValue)
c = (-minValue) / (maxValue - minValue)

from osgeo import gdal, gdal_array
import sys, os
file=os.getcwd() + '/doc.kml'
ds = gdal.Open(file)
if ds is None:
    print('Unable to open raster file')
    sys.exit(1)
band1 = ds.GetRasterBand(1)
band2 = ds.GetRasterBand(2)
band3 = ds.GetRasterBand(3)
rArray = band1.ReadAsArray()
gArray = band2.ReadAsArray()
bArray = band3.ReadAsArray()

data = np.zeros(rArray.shape)
i=0; j=0
for (rRow,gRow,bRow) in zip(rArray,gArray,bArray):
    j=0
    for (r,g,b) in zip(rRow,gRow,bRow):
        if r==35 and g==23 and b==27: data[i,j] = np.nan
        else:
            colorCode = symbols('x', real=True, positive=True)

            expr = 34.61 + colorCode * (1172.33 - colorCode * (10793.56 - colorCode * (33300.12 - colorCode * (38394.49 - colorCode * 14825.05))))
            sol_R = solve(expr - r, colorCode)
            #print("R: " + str(sol_R))

            expr = 23.31 + colorCode * (557.33 + colorCode * (1225.33 - colorCode * (3574.96 - colorCode * (1073.77 + colorCode * 707.56))))
            sol_G = solve(expr - g, colorCode)
            #print("G: " + str(sol_G))

            expr = 27.2 + colorCode * (3211.1 - colorCode * (15327.97 - colorCode * (27814 - colorCode * (22569.18 - colorCode * 6838.66))))
            sol_B = solve(expr - b, colorCode)
            #print("B: " + str(sol_B))

            '''iList = []
            arr = np.array(sol_R + sol_G + sol_B)
            arr2 = np.less_equal(np.abs(arr[:, np.newaxis] - arr[np.newaxis, :]), 0.1).astype(int)
            for iRow, iCol in np.transpose(np.triu_indices_from(arr2, k=1)):
                if arr2[iRow, iCol]:
                    iList.append([iRow, iCol])
                    print(iRow, iCol)
            ''' 

            colorCode = []
            delta=10
            for r in sol_R:
                for g in sol_G:
                    if np.absolute(r-g) < delta: 
                        colorCode.clear()
                        delta = np.absolute(r-g)
                        colorCode.append(r); colorCode.append(g)

            delta=10
            for b in sol_B:
                if np.absolute(b - Average(colorCode)) < delta: 
                    delta = np.absolute(b-Average(colorCode))
                    colorCodeValue = Average(colorCode+[b])

            value = symbols('x', real=True)
            expr = value*a + c
            sol = solve(expr - colorCodeValue, value)
            if sol[0] < minValue: sol[0] = minValue
            elif sol[0] > maxValue: sol[0] = maxValue
            #print(i,j)
            data[i,j] = sol[0]

        j+=1
    i+=1


# save array, using ds as a prototype
gdal_array.SaveArray(data.astype("float32"), os.getcwd() + '/power.tif', "GTIFF", ds)

ds = None

文件在这里:link

python-3.x numpy numpy-ndarray
1个回答
1
投票

[如果,就像您在评论中所说的那样,您只需要单个值(floatint或您拥有的值),一次遍历一个元素的数组的最佳方法是使用ndarray.ravel()ndarray.flatten()。两者之间的区别在于,ravel允许您在循环时修改原始数组,而flatten可以进行复制。对于您的情况,您不想修改数组,因此ravel is faster

请注意,无论如何都不会跟踪索引,因此代码中的ij毫无意义。

data = np.zeros(rArray.shape)
for r,g,b in zip(rArray.ravel(),gArray.ravel(),bArray.ravel()):
  res=2*r + g - b

编辑:由于您需要跟踪索引,因此最好像在索引中那样遍历两者。使用enumerate更具Python风格,这意味着您在编辑和调试代码时不必跟踪索引。我已经修改了循环块以使用这些变量以及变量,以便它们与循环变量不冲突。该代码对我有用。

for i,rows in enumerate(zip(rArray,gArray,bArray)):
    for j,(r,g,b) in enumerate(zip(*rows)):
        if r==35 and g==23 and b==27: data[i,j] = np.nan
        else:
            colorCode = symbols('x', real=True, positive=True)

            expr = 34.61 + colorCode * (1172.33 - colorCode * (10793.56 - colorCode * (33300.12 - colorCode * (38394.49 - colorCode * 14825.05))))
            sol_R = solve(expr - r, colorCode)
            #print("R: " + str(sol_R))

            expr = 23.31 + colorCode * (557.33 + colorCode * (1225.33 - colorCode * (3574.96 - colorCode * (1073.77 + colorCode * 707.56))))
            sol_G = solve(expr - g, colorCode)
            #print("G: " + str(sol_G))

            expr = 27.2 + colorCode * (3211.1 - colorCode * (15327.97 - colorCode * (27814 - colorCode * (22569.18 - colorCode * 6838.66))))
            sol_B = solve(expr - b, colorCode)
            #print("B: " + str(sol_B))

            colorCode = []
            delta=10
            for r_sol in sol_R:
                for g_sol in sol_G:
                    if np.absolute(r_sol-g_sol) < delta: 
                        colorCode.clear()
                        delta = np.absolute(r_sol-g_sol)
                        colorCode.append(r_sol); colorCode.append(g_sol)

            delta=10
            for b_sol in sol_B:
                if np.absolute(b_sol - Average(colorCode)) < delta: 
                    delta = np.absolute(b_sol-Average(colorCode))
                    colorCodeValue = Average(colorCode+[b_sol])

            value = symbols('x', real=True)
            expr = value*a + c 
            sol = solve(expr - colorCodeValue, value)
            if sol[0] < minValue:
                sol[0] = minValue
            elif sol[0] > maxValue:
                sol[0] = maxValue
            data[i,j] = sol[0]
© www.soinside.com 2019 - 2024. All rights reserved.