python openCVを用いたドロネー分割

Question

皆様の助けもあり画像から複数輪郭抽出を行うまでになりました．白の領域においてドロネー分割をおこないたいのですが，
どなたか解決できる方いらっしゃいましたら．よろしくお願いします．
python2.7,openCV4.0です．

#!/usr/bin/env python
#coding: utf-8

import csv,os,cv2,math,re
import numpy as np
import matplotlib.pyplot as plt
from matplotlib import pyplot
from operator import itemgetter
from dolfin import *
from mshr import *
from pylab import show,triplot

#読み込み,グレー設定
img = cv2.imread('/home/ubuntu/0116maps/map.pgm')
img_gray = cv2.cvtColor(img,cv2.COLOR_BGR2GRAY)

#smoothing
img_preprocessed = cv2.GaussianBlur(img_gray, (5, 5), 50)
img_preprocessed = cv2.GaussianBlur(img_preprocessed, (5, 5), 50)
img_preprocessed = cv2.GaussianBlur(img_preprocessed, (5, 5), 50)
img_preprocessed = cv2.GaussianBlur(img_preprocessed, (5, 5), 50)
img_preprocessed = cv2.GaussianBlur(img_preprocessed, (5, 5), 50)
img_preprocessed = cv2.GaussianBlur(img_preprocessed, (5, 5), 50)
img_preprocessed = cv2.GaussianBlur(img_preprocessed, (5, 5), 50)
img_preprocessed = cv2.GaussianBlur(img_preprocessed, (5, 5), 50)
img_preprocessed = cv2.GaussianBlur(img_preprocessed, (5, 5), 50)
img_preprocessed = cv2.GaussianBlur(img_preprocessed, (5, 5), 50)
img_preprocessed = cv2.GaussianBlur(img_preprocessed, (5, 5), 50)

#閾値処理,輪郭検索
_, white_binary = cv2.threshold(img_preprocessed, 220, 254, cv2.THRESH_BINARY)
white_contours, _ = cv2.findContours(white_binary, cv2.RETR_EXTERNAL, cv2.CHAIN_APPROX_SIMPLE)

#imageのコピー
white_and_contours = np.copy(img)

#最大輪郭値取得
min_white_area = 60
large_contours = [ cnt for cnt in white_contours if cv2.contourArea(cnt) > min_white_area ]

#外形の座標取得,最大値最小値取得
large_contours = np.array(large_contours)
large_contours_min = large_contours.min(axis = 1)
large_contours_max = large_contours.max(axis = 1)

#csv削除,書き込み
os.remove("/home/ubuntu/map_data/pixel_output.csv")
with open("/home/ubuntu/map_data/pixel_output.csv", 'a') as f:
    writer = csv.writer(f)
    writer.writerows(large_contours)

#pixel値の[]削除、,付け
os.remove("/home/ubuntu/map_data/pixel_surround.csv")
with open("/home/ubuntu/map_data/pixel_output.csv", 'rb') as f:
   reader = csv.reader(f)
   for row in reader:
       i = ' '.join(row)
       i = i.replace(' ',',')
       i = i.replace('],[','')
       i = i.replace('[','')
       i = i.replace(']','\n')
       f = open("/home/ubuntu/map_data/pixel_surround.csv", 'a')
       f.write(i)
       f.close()

os.remove("/home/ubuntu/map_data/minmaxoutput.csv")
with open("/home/ubuntu/map_data/minmaxoutput.csv", 'a') as f:
    writer = csv.writer(f)
    writer.writerows(large_contours_min)
    writer.writerows(large_contours_max)

os.remove("/home/ubuntu/map_data/minmax_coordinate.csv")
with open("/home/ubuntu/map_data/minmaxoutput.csv", 'rb') as f:
    reader = csv.reader(f)
    for row in reader:
      i = ' '.join(row)
      i = i.replace('[','')
      i = i.replace(']',',')
      i = i.replace(' ',',')
      f = open("/home/ubuntu/map_data/minmax_coordinate.csv", 'a')
      f.write(i)
      f.close()

#map値読み取りやすく変換
os.remove("/home/ubuntu/map_data/map_value.csv")
with open("/home/ubuntu/0116maps/map.yaml", 'rb') as f:
  reader = csv.reader(f)
  for words in reader:
     i = ' '.join(words)
     match = re.search(r': ', i)
     if match:
        i = i[match.end():]
	i = i.strip('[]')
	i = i.replace('  ',',')
	f = open("/home/ubuntu/map_data/map_value.csv", 'a')
	f.write(i + ",")
	f.close()

#格子分割準備
with open("/home/ubuntu/map_data/pixel_surround.csv", 'rb') as f:
   reader = csv.reader(f)
   count = 0
   for pixel in reader:
      count+=1

j = count-1
k = 0
poin = range(j)
domain_vertices = []

with open("/home/ubuntu/map_data/pixel_surround.csv", 'rb') as f:
   reader = csv.reader(f)
   for pix in reader:
     if k<j:
        poin[k] = [Point(int(pix[0]),int(pix[1]))]
        k+=1

for m in reversed(range(j)):
        domain_vertices += poin[m]

#格子分割
domain = Polygon(domain_vertices)
mesh = generate_mesh(domain,0.1)
coords = mesh.coordinates()

os.remove("/home/ubuntu/map_data/coords.csv")
with open("/home/ubuntu/map_data/coords.csv", 'a') as f:
    writer = csv.writer(f)
    writer.writerows(coords)
i
#メッシュ領域plot
triplot(coords[:,0], coords[:,1], triangles=mesh.cells())
plt.imshow(white_and_contours)
plt.show()

Comments

[tiitoi]

ソースコードも載せていただけますか？
コードがないと、どう変えればいいのか回答できないです。

[gonzoshira]

すいません．添付しました！

[tiitoi]

mshr は使ったことがないので、自分は回答できないです。すいません。

[gonzoshira]

私もこのソースコードを頂いたばかりでして．．．
質問なのですが，mshrとdolfinとはなんあのでしょうか？

[tiitoi]

有限要素法の領域分割に使うライブラリのようです。
使ったことがないので、私も詳しいことはわかりません。
https://people.sc.fsu.edu/~jburkardt/py_src/dolfin/dolfin.html
その頂いたコードを作成した人に聞くのがいいかもしれません。

[gonzoshira]

いろいろご指導いただきありがとうござます

Answer 1

以下、直接の回答ではありません。参考まで。

Triangleというライブラリを見つけました。
入出力データ形式がちょっとよく分かりませんが、比較的簡単に「穴」を除外して分割できるようです。
使い方や分割イメージはExamplesやdataを見ると何となく分かるかと思います。

import numpy as np
import matplotlib.pyplot as plt
import triangle as tr

#box = tr.get_data('box') # サンプルデータ、穴の開いた四角形の読込

# 独自に与えるデータ。最低限以下があればよいようだ。
box = {'vertices': np.array([
    [0.5, 0.5], [0., 3.],[3. , 0. ],[3. , 3. ],
    [1. , 1. ],[1. , 2. ],[2. , 1. ],[2. , 2. ]]),
 'segments': np.array([[0, 1],[4, 6],[6, 7],[7, 5],[5, 4]]),
 'holes': np.array([[1.5, 1.5]])
}

print(box)

# 分割結果は主に'triangles'に格納されるようだ
t = tr.triangulate(box, 'pc')
#t = tr.triangulate(box, 'pcqa0.1') # より細かく
print(t)

tr.compare(plt, box, t) # 比較描画
plt.savefig('ret.png')

Comments

[gonzoshira]

貴重な情報ありがとうございます！！