翻譯

2014-07-12 171 views
1

我有我生成使用翻譯

A1, M1, A2, M2, A3, M3, E = meshgrid(Grid.aGrid, Grid.mGrid, Grid.aGrid, Grid.mGrid, Grid.aGrid, Grid.mGrid, Grid.eGrid, indexing='ij') 

,其中Grid.aGrid使用linspace(aMin, aMax, nA)產生形狀(A1, M1, A2, M2, A3, M3, E)的網格,並且類似地對於其它的網格。

考慮一些Z = f(A1, ...),其中f()將標記一些網格點爲無關緊要。爲簡單起見,假設它是

Z = A1 + A2 + A3 
Z[Z < 0] = NaN 

考慮Z[0, 1, 2, 3, 4, 5, 6]。它包含對應於實際值(aGrid[0], mGrid[1], aGrid[2], mGrid[3], aGrid[4], mGrid[5], eGrid[6])的值。這正是我嘗試實現未標記的f()Z所有點:

我想創建一個字典

foo = {z1, z2, z3, ... zn} 

其中z1等都是那種

z1 = (aGrid[0], mGrid[1], aGrid[2], mGrid[3], aGrid[4], mGrid[5], eGrid[6]) 

,這是Z內部對應z1位置的網格值。

我想出了:

aGrid = arange(0, 10) 
mGrid = arange(100, 110) 
eGrid = arange(1000, 1200) 

A,M,E = meshgrid(aGrid, mGrid, eGrid, indexing='ij') 

# contains the grid index 
Z = (A + M + E).astype(float) 
Z[A < 3] = nan 
# will contain the actual values, as tuples 
Z2 = {} 

for i, idx in enumerate(ndindex(Z.shape)): 
    a = aGrid[idx[0]] 
    m = mGrid[idx[1]] 
    e = eGrid[idx[2]] 

    if isnan(Z[idx]): 
     Z2[i] = NaN 
    else: 
     Z2[i] = (a, m, e) 

效率是關鍵。有沒有更快更清潔的方式可以實現這一目標?任何使用字典的替代方法?

我特別不喜歡我必須寫下aGrid[idx[0]]等。是否有可能讓算法更一般?沿的

for i, idx in enumerate(ndindex(Z.shape)): 
    # some magic happens here. What exactly? 
    someMagicList = magic(aGrid, mGrid, eGrid) 
    Z2[i] = someMagicList[idx] 

回答

1

broadcast_arrays()使用的線,並且將結果Z2一些事情是具有形狀(20000, 3)陣列。

import numpy as np 

aGrid = np.arange(0, 10, dtype=float) 
mGrid = np.arange(100, 110, dtype=float) 
eGrid = np.arange(1000, 1200, dtype=float) 

A,M,E = np.meshgrid(aGrid, mGrid, eGrid, indexing='ij') 

# contains the grid index 
Z = (A + M + E).astype(float) 
Z[A < 3] = np.nan 

grids = [A, M, E] 
grid_bc = np.broadcast_arrays(*grids) 
Z2 = np.column_stack([g.ravel() for g in grid_bc]) 
Z2[np.isnan(Z.ravel())] = np.nan 

print Z2[5900], Z2[6000] 
+0

這看起來比我所做的更有效率的提高。正如我正在試圖將其寫入我的腳本中;如果可能的話,我會做一個簡短的後續行動(我認爲這對於一個新問題來說太短了):我如何(有效地)獲取那些不是「NaN」的子集? 'stateRelevant = [z for z in states如果不是np.isnan(z.any())]'把它們全部放進去了。顯然,'isnan()'沒有用來選擇那些具有'(nan,.. )在他們身上。 – FooBar

+0

我不明白你的問題,'州'是什麼?你是指'Z [〜np.isnan(Z)]'還是'Z2 [〜np.isnan(Z)]'? – HYRY

+0

對不起,我的意思是'Z2'。使用'Z'作爲'isnan()'的指示器很好,但如果我只剩下'Z2',該怎麼辦?我打算只通過處理它的函數傳遞'Z2'。 – FooBar