かっこつけて英語のタイトルにしようとして果てる、の巻。
英語のタイトル付けることが多いのは、元ネタが StackOverflow のことが多いからなの。そいで、StackOverflow からのネタが多いのは、これとかこれが理由。
つーか、日本語でも伝わりにくいなぁ…。立体地形図を作るために色んな領域を扱ってるとね、「海なんかいらないのに!」となることが多いのね。陸域だけに興味があるわけだからさ。例えば伊豆半島南端ぴったりを南端に、銚子をぴったり東端に、とするのに、地理院地図WEBや Google Map で調べるのは、なかなかに苦痛だ。
海データは NaN にしている、ワタシは。つまり、入力のこれ：
1 >>> import numpy as np
2 >>> a = np.full((10, 11), np.nan)
3 >>> a = np.triu(a) + np.arange(10*11).reshape((10, 11))
4 >>> a[-2:,:2] = np.nan
5 >>> a[-4:-2,2:4] = np.nan
6 >>> a
7 array([[ nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan],
8 [ 11., nan, nan, nan, nan, nan, nan, nan, nan, nan, nan],
9 [ 22., 23., nan, nan, nan, nan, nan, nan, nan, nan, nan],
10 [ 33., 34., 35., nan, nan, nan, nan, nan, nan, nan, nan],
11 [ 44., 45., 46., 47., nan, nan, nan, nan, nan, nan, nan],
12 [ 55., 56., 57., 58., 59., nan, nan, nan, nan, nan, nan],
13 [ 66., 67., nan, nan, 70., 71., nan, nan, nan, nan, nan],
14 [ 77., 78., nan, nan, 81., 82., 83., nan, nan, nan, nan],
15 [ nan, nan, 90., 91., 92., 93., 94., 95., nan, nan, nan],
16 [ nan, nan, 101., 102., 103., 104., 105., 106., 107., nan, nan]])
から、これを抽出したいわけである：
1 >>> # 何か操作…
2 >>> a
3 array([[ 11., nan, nan, nan, nan, nan, nan, nan, nan],
4 [ 22., 23., nan, nan, nan, nan, nan, nan, nan],
5 [ 33., 34., 35., nan, nan, nan, nan, nan, nan],
6 [ 44., 45., 46., 47., nan, nan, nan, nan, nan],
7 [ 55., 56., 57., 58., 59., nan, nan, nan, nan],
8 [ 66., 67., nan, nan, 70., 71., nan, nan, nan],
9 [ 77., 78., nan, nan, 81., 82., 83., nan, nan],
10 [ nan, nan, 90., 91., 92., 93., 94., 95., nan],
11 [ nan, nan, 101., 102., 103., 104., 105., 106., 107.]])
ググってもピンポイントではヒットせんなぁ。
仕方ないので自力でみたが。もうね、一撃でとかそんなん無理みたいで。まっとうにやるしかなかったわ：
1 >>> def get_bound(arr2d, f=np.isnan):
2 ... #
3 ... def _sl(arr2d, i, axis):
4 ... if axis == 0:
5 ... return arr2d[i,:]
6 ... else:
7 ... return arr2d[:,i]
8 ... #
9 ... def _gb(arr2d, axis, f):
10 ... idx_l = 0
11 ... for i in range(a.shape[axis]):
12 ... idx_l = i
13 ... if not np.all(f(_sl(a, i, axis))):
14 ... break
15 ... idx_h = a.shape[axis]
16 ... for i in range(a.shape[axis], 1, -1):
17 ... idx_h = i
18 ... if not np.all(f(_sl(a, i - 1, axis))):
19 ... break
20 ... return slice(idx_l, idx_h)
21 ... #
22 ... return _gb(arr2d, 0, f), _gb(arr2d, 1, f)
23 ...
24 >>> get_bound(a)
25 (slice(1, 10, None), slice(0, 9, None))
26 >>> a[get_bound(a)]
27 array([[ 11., nan, nan, nan, nan, nan, nan, nan, nan],
28 [ 22., 23., nan, nan, nan, nan, nan, nan, nan],
29 [ 33., 34., 35., nan, nan, nan, nan, nan, nan],
30 [ 44., 45., 46., 47., nan, nan, nan, nan, nan],
31 [ 55., 56., 57., 58., 59., nan, nan, nan, nan],
32 [ 66., 67., nan, nan, 70., 71., nan, nan, nan],
33 [ 77., 78., nan, nan, 81., 82., 83., nan, nan],
34 [ nan, nan, 90., 91., 92., 93., 94., 95., nan],
35 [ nan, nan, 101., 102., 103., 104., 105., 106., 107.]])
36 >>>
37 >>> # ややエレガントにみえるバージョン
38 >>> # (ndimsが1でもいける(確認済み)。3でもいけるはず? (3は未確認。))
39 >>> def get_bound(arr, f=np.isnan):
40 ... #
41 ... def _gb(arr, axis, f):
42 ... _sl = [slice(None), ] * arr.ndim
43 ... idx_l = 0
44 ... for i in range(arr.shape[axis]):
45 ... idx_l = i
46 ... _sl[axis] = slice(i, i + 1)
47 ... if not np.all(f(arr[_sl])):
48 ... break
49 ... idx_h = arr.shape[axis]
50 ... for i in range(arr.shape[axis], 1, -1):
51 ... idx_h = i
52 ... _sl[axis] = slice(i - 1, i)
53 ... if not np.all(f(arr[_sl])):
54 ... break
55 ... return slice(idx_l, idx_h)
56 ... #
57 ... return [_gb(arr, axis, f) for axis in range(arr.ndim)]
58 ...
59 >>> get_bound(a)
60 [slice(1, 10, None), slice(0, 9, None)]
61 >>> a[get_bound(a)]
62 array([[ 11., nan, nan, nan, nan, nan, nan, nan, nan],
63 [ 22., 23., nan, nan, nan, nan, nan, nan, nan],
64 [ 33., 34., 35., nan, nan, nan, nan, nan, nan],
65 [ 44., 45., 46., 47., nan, nan, nan, nan, nan],
66 [ 55., 56., 57., 58., 59., nan, nan, nan, nan],
67 [ 66., 67., nan, nan, 70., 71., nan, nan, nan],
68 [ 77., 78., nan, nan, 81., 82., 83., nan, nan],
69 [ nan, nan, 90., 91., 92., 93., 94., 95., nan],
70 [ nan, nan, 101., 102., 103., 104., 105., 106., 107.]])
こんなん、かなり一般的なタスクなんじゃないかって気がするんだけどなぁ…。例えば画像処理でまわりの黒い部分を取り除きたい、とかさ。
なお、「自力た」言うても、これはヒントにした。
ほいだば、実践だす。
ラフに選択した範囲ではこうなるとして：
1 # elevs はどこかからもってきた標高データ (上の絵)
2 elevs = elevs[get_bound(elevs)] # 海ばっかのとこはいらない
うむ。よか。
ちょっとだけバリエーション。elevs.T (転置) と np.fliplr, np.flipud だけで欠けの方向のバリエーションをすぐに作れるのでね。
↓
↓
↓
↓
↓
↓
↓
↓
↓
島だろうとヘコたれない。
↓
↓
↓
