From e421cea33122a735ad4f7fb0bacfcca12a4e00b0 Mon Sep 17 00:00:00 2001 From: Evgeni Burovski Date: Tue, 12 Jan 2016 16:20:32 +0000 Subject: [PATCH] ENH: start adding a python version of Skeletonize3D code --- skimage/morphology/data/3/pngs/1.png | Bin 0 -> 251 bytes skimage/morphology/data/3/pngs/2.png | Bin 0 -> 248 bytes skimage/morphology/data/3/pngs/3.png | Bin 0 -> 184 bytes skimage/morphology/data/3/pngs/5.png | Bin 0 -> 178 bytes skimage/morphology/data/3/pngs/6.png | Bin 0 -> 178 bytes skimage/morphology/data/3/pngs/7.png | Bin 0 -> 176 bytes skimage/morphology/data/3/stack.tif | Bin 0 -> 15187 bytes skimage/morphology/data/3/stack_fiji.tif | Bin 0 -> 15187 bytes skimage/morphology/data/4/1_fiji.png | Bin 0 -> 94 bytes skimage/morphology/data/4/pngs/1.png | Bin 0 -> 149 bytes skimage/morphology/data/4/pngs/2.png | Bin 0 -> 137 bytes skimage/morphology/data/4/pngs/3.png | Bin 0 -> 137 bytes skimage/morphology/data/4/pngs/4.png | Bin 0 -> 137 bytes skimage/morphology/data/4/pngs/5.png | Bin 0 -> 137 bytes skimage/morphology/data/4/pngs/6.png | Bin 0 -> 137 bytes skimage/morphology/data/4/stack.tif | Bin 0 -> 1795 bytes skimage/morphology/data/4/stack_fiji.tif | Bin 0 -> 1795 bytes skimage/morphology/data/cross.png | Bin 0 -> 149 bytes skimage/morphology/data/cross.txt | 24 + skimage/morphology/data/cross_fiji.txt | 24 + skimage/morphology/data/loop.png | Bin 0 -> 184 bytes skimage/morphology/data/loop.txt | 24 + skimage/morphology/data/loop_fiji.txt | 24 + skimage/morphology/data/strip.png | Bin 0 -> 131 bytes skimage/morphology/data/strip.txt | 12 + skimage/morphology/data/strip_fiji.txt | 12 + skimage/morphology/data/two-hole.png | Bin 0 -> 176 bytes skimage/morphology/data/two-hole.txt | 16 + skimage/morphology/data/two-hole_fiji.txt | 16 + skimage/morphology/skel.py | 678 ++++++++++++++++++++++ skimage/morphology/test_simple_2d.py | 84 +++ skimage/morphology/trivial.py | 73 +++ 32 files changed, 987 insertions(+) create mode 100644 skimage/morphology/data/3/pngs/1.png create mode 100644 skimage/morphology/data/3/pngs/2.png create mode 100644 skimage/morphology/data/3/pngs/3.png create mode 100644 skimage/morphology/data/3/pngs/5.png create mode 100644 skimage/morphology/data/3/pngs/6.png create mode 100644 skimage/morphology/data/3/pngs/7.png create mode 100644 skimage/morphology/data/3/stack.tif create mode 100644 skimage/morphology/data/3/stack_fiji.tif create mode 100644 skimage/morphology/data/4/1_fiji.png create mode 100644 skimage/morphology/data/4/pngs/1.png create mode 100644 skimage/morphology/data/4/pngs/2.png create mode 100644 skimage/morphology/data/4/pngs/3.png create mode 100644 skimage/morphology/data/4/pngs/4.png create mode 100644 skimage/morphology/data/4/pngs/5.png create mode 100644 skimage/morphology/data/4/pngs/6.png create mode 100644 skimage/morphology/data/4/stack.tif create mode 100644 skimage/morphology/data/4/stack_fiji.tif create mode 100644 skimage/morphology/data/cross.png create mode 100644 skimage/morphology/data/cross.txt create mode 100644 skimage/morphology/data/cross_fiji.txt create mode 100644 skimage/morphology/data/loop.png create mode 100644 skimage/morphology/data/loop.txt create mode 100644 skimage/morphology/data/loop_fiji.txt create mode 100644 skimage/morphology/data/strip.png create mode 100644 skimage/morphology/data/strip.txt create mode 100644 skimage/morphology/data/strip_fiji.txt create mode 100644 skimage/morphology/data/two-hole.png create mode 100644 skimage/morphology/data/two-hole.txt create mode 100644 skimage/morphology/data/two-hole_fiji.txt create mode 100644 skimage/morphology/skel.py create mode 100644 skimage/morphology/test_simple_2d.py create mode 100644 skimage/morphology/trivial.py diff --git a/skimage/morphology/data/3/pngs/1.png b/skimage/morphology/data/3/pngs/1.png new file mode 100644 index 0000000000000000000000000000000000000000..da940f11a60ff0c8e90ff61a0ffe6d000ad573b8 GIT binary patch literal 251 zcmeAS@N?(olHy`uVBq!ia0vp^Wz#oq6^RZsftubPXO7TIfX ztPA1UWge4v-0y4ChsE0$%~+UycHNec%re7EU*!JC@n^f{dnamtv)-}P)7C%q{dB+L t&vTmHo|WlNlYC~ua5+EwNq+q`_CJqpTPD6SUI}y;gQu&X%Q~loCIGAdUiSb1 literal 0 HcmV?d00001 diff --git a/skimage/morphology/data/3/pngs/2.png b/skimage/morphology/data/3/pngs/2.png new file mode 100644 index 0000000000000000000000000000000000000000..ada13a9ee12f81c9f75611f195dc52270cc6b870 GIT binary patch literal 248 zcmeAS@N?(olHy`uVBq!ia0vp^WU`533hbv~qyzDus+}e6r`e^mMz<_;bE@$^Jho6WudFkV@fBMZ= qmFxaV&D$8h#<}GGMyvn7XY)7c+0A7P650TC6oaR$pUXO@geCy`E?md} literal 0 HcmV?d00001 diff --git a/skimage/morphology/data/3/pngs/3.png b/skimage/morphology/data/3/pngs/3.png new file mode 100644 index 0000000000000000000000000000000000000000..00ee7d1f35fd2a8a5d415ba1be23eaf78328f066 GIT binary patch literal 184 zcmeAS@N?(olHy`uVBq!ia0vp^WviwG%p$wj`elF{r5}E)ar#q|w literal 0 HcmV?d00001 diff --git a/skimage/morphology/data/3/pngs/5.png b/skimage/morphology/data/3/pngs/5.png new file mode 100644 index 0000000000000000000000000000000000000000..8378036e08c6dbb5b5c81e1cea2d8ef961825c35 GIT binary patch literal 178 zcmeAS@N?(olHy`uVBq!ia0vp^W;vq2^c)t; zo-ml*xR5U|i>YPW{{wXlE}NZSSfz2~9{gTe~DWM4f6qGbC literal 0 HcmV?d00001 diff --git a/skimage/morphology/data/3/pngs/6.png b/skimage/morphology/data/3/pngs/6.png new file mode 100644 index 0000000000000000000000000000000000000000..3c6051aedc05befee077d00896d023f62c02f212 GIT binary patch literal 178 zcmeAS@N?(olHy`uVBq!ia0vp^W;vq2^c)t; zo-n#QPqFbQiHHKJLoIHLyB|Li=AOZM<1`P0zAv}Z Unsf7B0gYquboFyt=akR{0A9>Cy8r+H literal 0 HcmV?d00001 diff --git a/skimage/morphology/data/3/pngs/7.png b/skimage/morphology/data/3/pngs/7.png new file mode 100644 index 0000000000000000000000000000000000000000..61eb70015833d281f406b02d04e2fd2c0a720c23 GIT binary patch literal 176 zcmeAS@N?(olHy`uVBq!ia0vp^W;vq2^c)t; zo-n#QPqFbQiHHKJLp^SbyB|Li=H__C%wX_^ThTyN RrT}OdgQu&X%Q~loCIHSVGu!|G literal 0 HcmV?d00001 diff --git a/skimage/morphology/data/3/stack.tif b/skimage/morphology/data/3/stack.tif new file mode 100644 index 0000000000000000000000000000000000000000..28752135debf71383648d205de0b17e81f45646a GIT binary patch literal 15187 zcmeI3%}x|S5XXB^5jV!f81>}UgCdFs{7AUNCVTea2@MjGRY-`BqTYNHuOvKx2M->7 z4Ym5Kx~iwAr)B_yu}F6K*T=7_yJmJ~A=%l?+%v{3nXBfPxd^I(mpw$i-ZZF}n+0>h zg_j^_F=$9`u9*e!--p+i=6cL;9yB*e&aI%AA^&#Fxw3im$?40$9v^*n{Wn4XH0Ifl zt#{kIJJXHTm9!_Zhf&%(;o+j|y{tQ@<%g2KR#N;IUv6M?xFGg~S?7>K; zpG0s}u?UeQ0p?N&Vx`A2Cp8PL`mqRQ=Qg-$d8!k|xhh5>(4<}s^Rck|5=QEv3{V`J zH8JK(B{M0c2$w~ZiHcCjVNld8Vgf3wC)-!Me(IDJBmt0aYupvA?ntFs4eG=kM9yBe zoOalI>b?BVLRA7-nKfW$sG1QXs=!oR#er~+t)^)3{Q7oCfrJrfzqD9?mS*S#bOI-u zfO~%qsF2>z(zhrlg#8MT;D%77`tBnPvXcmo>QEv?k_4DbA&8Y8%be6KwCcwqlnpJo zX?dy>#knd*A<(2=4fC-$>jFmVpbStPnl&-zOC>WYqzIQqlZlE@$YD^_EMfvGt0&u6 zyMF4F6(j+WZfo2XtnNsqSq@cjCANP&bAXTP*qf0kzG1kUpW+{Hg|f#D4^kYYEU0qmF+*6!c6*5Ws`IFWT5 zb;B5So3TSpcDB2CNSpC=@|Mfl8wz{vRsk{7_p>sBqc4cc_*B}sEeTmi7pL7 zc?m~cz4M?&1fklWqdF=6#{#QS{cO(*f`u#L5p?&6ze)8hMZ-D)oq$e2C(xe2nA@~t z$OEEXwdnC!F;sH`RttB;Mp7Ns)P3p`$gT!KA`eVoxQ2IHlG*$KK1at!b(gY@_|E0#s2g#Qjjd}hG pzK3|~U#YzDRoMOeS#-a1N!Ig%pX-K&JzFt literal 0 HcmV?d00001 diff --git a/skimage/morphology/data/3/stack_fiji.tif b/skimage/morphology/data/3/stack_fiji.tif new file mode 100644 index 0000000000000000000000000000000000000000..9373767547c7d86a4192e6072677d61b4daaa679 GIT binary patch literal 15187 zcmeI3%SyvQ6oyagJ%}KpuH3rtQg5_gpnF^B+J&nswG^TiijSi1d=s~V58%RuOJBo! zk~U*5lXhloLR45TdPj}Ry zLlGT`s34-Eh>k_X4RC;f7J&{|yQk{WagT1X|9g-$TBl?XBSkx6PKsnE*Yfl5tGq?IsqFq8FB^_eoHkk~PC!ULByo0RX; zPMA$di#sHpmnMb84dNBj=Q*T_s>HSvsfrDKv%ifhRw7BTN>9-3;hRZ1A_W2+6gko~VPY9BL|YIoFu00@8p2uwHu ze2h%^4Ib6yMIOA+Ffm|s00ck)1VCW){}hIMG~NUf|EDGLT`4jg??vP4K>!3m z00ck)1VA7|0{9-4A#7NI00<11fIa+|8tn&_ml}yqFZn&hp}$hO@~Y_mnTGqujGFt| UMH=p#RYmuYG~D+d)!d1G0j7(b{r~^~ literal 0 HcmV?d00001 diff --git a/skimage/morphology/data/4/1_fiji.png b/skimage/morphology/data/4/1_fiji.png new file mode 100644 index 0000000000000000000000000000000000000000..03a7b79a4570ba4c6a89d04c3bc351c88bfcb1e3 GIT binary patch literal 94 zcmeAS@N?(olHy`uVBq!ia0vp^JRr;gBp8b2n5}`7s;7%%NX4z>gaiEllnwtYSNxy& qBw)pVo*6bRH|p6^l>hxdo^Q^ewny;2b@lEEAdQ}`elF{r5}E)*&K%YN literal 0 HcmV?d00001 diff --git a/skimage/morphology/data/4/pngs/1.png b/skimage/morphology/data/4/pngs/1.png new file mode 100644 index 0000000000000000000000000000000000000000..d03660d2d80e785af86148e963b0466643aab2a4 GIT binary patch literal 149 zcmeAS@N?(olHy`uVBq!ia0vp^JRr=-3?y@QpIZW?n2Vh}LpV4%Za?&Y0OYa-_=LCu z>Hq)#FTFYcE|A4p;1OBOz`%C|gc+x5^GO2*)jVAsLpWw8CnO{%BqmHqPS}v2;9+z4 o!vjW+sWu54><=i^H!$)rOjh7HzImO`Y@iMXPgg&ebxsLQ0NPY4f&c&j literal 0 HcmV?d00001 diff --git a/skimage/morphology/data/4/pngs/2.png b/skimage/morphology/data/4/pngs/2.png new file mode 100644 index 0000000000000000000000000000000000000000..efa5c10492ba48e3074e296370deb5054d498853 GIT binary patch literal 137 zcmeAS@N?(olHy`uVBq!ia0vp^JRr=-3?y@QpIZW?n2Vh}LpV4%Za?&Y0OYa-_=LCu z>Hq)#FTFYcE|A4p;1OBOz`%C|gc+x5^GO2*r9536LpWw8CnTgbG_#1cvGefo?~oLd bn83)Oe~D?bz!bqIpehDWS3j3^P6Hq)#FTFYcE|A4p;1OBOz`%C|gc+x5^GO2*r9536LpWw8CnO{{FthWuG4t^7OGpYy bOkiZtn!psGcy_^QpehDWS3j3^P6Hq)#FTFYcE|A4p;1OBOz`%C|gc+x5^GO2*r9536LpWw8CnO{{FthWuG4t^7OGpYy bOkiZtn!psGcy_^QpehDWS3j3^P6Hq)#FTFYcE|A4p;1OBOz`%C|gc+x5^GO2*r9536LpWw8CnO{{FthWuG4t^7OGpYy bOkiZtn!psGcy_^QpehDWS3j3^P6Hq)#FTFYcE|A4p;1OBOz`%C|gc+x5^GO2*r9536LpWw8CnO{{FthWuG4t^7OGpYy bOkiZtn!psGcy_^QpehDWS3j3^P6%#Pp_QBq#U3e8-pQ(9zy{` z9z!~V5u9ZVXPLlRrf`-SoCP!h1^kEe{sRG&17ZRJhzBEq9H5Yu25LT#nU3X(J<`3=cw@UR+n4h3&GWil|dnZk7ze}9$Z0xkKojI3C+goX3`(NU0XJ7yT2!@p& literal 0 HcmV?d00001 diff --git a/skimage/morphology/data/4/stack_fiji.tif b/skimage/morphology/data/4/stack_fiji.tif new file mode 100644 index 0000000000000000000000000000000000000000..c0d182fd31f51f6b1ae14876343a46bd4cf834dd GIT binary patch literal 1795 zcmebEWzb?^VBlcjVfe?u0%S1)F^GV%c^Day*h~z}U}+AZ8a5~!sF;zDfe9$C4aD;p z1!3wLbs0sG)QCaZJV5p0Ff|Mlj!grz*@1db9-9uf-vh{g3UorRXKrG8s+X;yo~c19 zS0;#4Y-`3%#Pp_QBq#U3e8-pQ(9zy{` z9z!~V5u9ZVXPLlRrf`-SoCP!hC-@JS`Tw5*#sSkX4p;#u2O}(d(|4_WO0;q$*)78&qol`;+03>EBBLDyZ literal 0 HcmV?d00001 diff --git a/skimage/morphology/data/cross.txt b/skimage/morphology/data/cross.txt new file mode 100644 index 00000000..2118f7f1 --- /dev/null +++ b/skimage/morphology/data/cross.txt @@ -0,0 +1,24 @@ +0 0 0 0 0 255 255 255 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 255 255 255 255 255 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 255 255 255 255 255 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 255 255 255 255 255 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 255 255 255 255 255 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 255 255 255 255 255 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 255 255 255 255 255 255 255 255 255 255 255 255 255 255 255 255 255 255 255 0 0 0 +0 255 255 255 255 255 255 255 255 255 255 255 255 255 255 255 255 255 255 255 255 0 0 0 +0 255 255 255 255 255 255 255 255 255 255 255 255 255 255 255 255 255 255 255 255 255 0 0 +0 255 255 255 255 255 255 255 255 255 255 255 255 255 255 255 255 255 255 255 255 255 0 0 +0 0 255 255 255 255 255 255 255 255 255 255 255 255 255 255 255 255 255 255 255 255 0 0 +0 0 0 0 0 255 255 255 255 255 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 255 255 255 255 255 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 255 255 255 255 255 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 255 255 255 255 255 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 255 255 255 255 255 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 255 255 255 255 255 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 255 255 255 255 255 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 255 255 255 255 255 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 255 255 255 255 255 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 255 255 255 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 diff --git a/skimage/morphology/data/cross_fiji.txt b/skimage/morphology/data/cross_fiji.txt new file mode 100644 index 00000000..c4ea69a9 --- /dev/null +++ b/skimage/morphology/data/cross_fiji.txt @@ -0,0 +1,24 @@ +0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 255 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 255 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 255 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 255 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 255 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 255 255 255 255 0 255 255 255 255 255 255 255 255 255 255 255 255 0 0 0 0 +0 0 0 0 0 0 0 255 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 255 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 255 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 255 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 255 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 255 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 255 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 255 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 255 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 diff --git a/skimage/morphology/data/loop.png b/skimage/morphology/data/loop.png new file mode 100644 index 0000000000000000000000000000000000000000..05fdc25955801664638367107253540ee3885375 GIT binary patch literal 184 zcmeAS@N?(olHy`uVBq!ia0vp^l0YoM$P6S+tB1s;*b3=DinK$vl=HlH+5(9_e!F@$4gaza9aMw$crH~xkr^#UKy zGkDZ62=(wcuzYUhcX;%mS&Hr0V_t{OYK5cq4KqIS_Z@W!H)#L> literal 0 HcmV?d00001 diff --git a/skimage/morphology/data/loop.txt b/skimage/morphology/data/loop.txt new file mode 100644 index 00000000..277e1b47 --- /dev/null +++ b/skimage/morphology/data/loop.txt @@ -0,0 +1,24 @@ +0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 0 0 255 255 255 255 255 255 255 255 255 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 255 255 255 255 255 255 255 255 255 255 255 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 255 255 255 255 255 255 255 0 0 0 255 255 0 0 0 0 0 +0 0 0 0 0 0 0 0 255 255 255 255 0 0 0 0 0 0 255 255 0 0 0 0 0 +0 0 0 0 0 0 0 255 255 255 255 0 0 0 0 0 0 0 255 255 0 0 0 0 0 +0 0 0 0 0 0 255 255 255 255 0 0 0 0 0 0 255 255 255 255 0 0 0 0 0 +0 0 0 0 0 0 255 255 255 255 0 0 0 0 0 255 255 255 255 0 0 0 0 0 0 +0 0 0 0 0 0 0 255 255 255 0 0 0 255 255 255 255 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 255 255 0 255 255 255 255 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 255 255 255 255 255 255 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 255 255 255 255 255 255 255 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 255 255 255 255 255 255 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 255 255 255 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 255 255 255 255 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 0 255 255 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 0 255 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 0 255 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 0 255 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 diff --git a/skimage/morphology/data/loop_fiji.txt b/skimage/morphology/data/loop_fiji.txt new file mode 100644 index 00000000..4d38bfac --- /dev/null +++ b/skimage/morphology/data/loop_fiji.txt @@ -0,0 +1,24 @@ +0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 0 0 0 255 255 255 255 255 255 255 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 0 0 255 0 0 0 0 0 0 0 255 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 0 255 0 0 0 0 0 0 0 0 255 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 255 0 0 0 0 0 0 0 0 0 255 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 255 0 0 0 0 0 0 0 0 0 255 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 255 0 0 0 0 0 0 0 255 255 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 255 0 0 0 0 0 0 255 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 255 0 0 0 0 255 255 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 255 0 0 0 255 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 0 255 255 255 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 0 255 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 0 255 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 0 255 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 0 255 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 0 255 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 0 255 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 0 255 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 diff --git a/skimage/morphology/data/strip.png b/skimage/morphology/data/strip.png new file mode 100644 index 0000000000000000000000000000000000000000..74e9b3e49c0730a12cbf27c91e9dd825acd61822 GIT binary patch literal 131 zcmeAS@N?(olHy`uVBq!ia0vp^+(693$P6U=(y!hHQq09po*^6@9Je3(KLB#s0(?ST z|Ns9FWQHEPTnD5$3p^r=85sDEfH31!Z9ZwBpopi7V+hC0|22Vm7w!=Se(5=Wj}+%v){v!G=Hi!Tzopr0FyU3r~m)} literal 0 HcmV?d00001 diff --git a/skimage/morphology/data/two-hole.txt b/skimage/morphology/data/two-hole.txt new file mode 100644 index 00000000..680e7beb --- /dev/null +++ b/skimage/morphology/data/two-hole.txt @@ -0,0 +1,16 @@ +0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 0 255 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 255 255 255 255 255 0 0 0 0 0 +0 0 0 0 0 255 255 255 255 255 255 255 255 255 0 255 0 +0 0 0 255 255 255 255 255 255 255 255 255 255 255 255 255 0 +0 0 255 255 255 255 0 0 255 255 255 0 0 255 255 255 0 +0 0 255 255 255 0 0 0 255 255 255 0 0 255 255 255 0 +0 0 0 255 255 0 0 0 255 255 255 0 0 255 255 255 0 +0 0 0 255 255 255 0 0 255 255 255 0 0 255 255 255 0 +0 0 0 0 255 255 0 0 255 255 255 0 0 255 255 255 0 +0 0 0 0 255 255 255 255 255 255 255 255 255 255 255 255 0 +0 0 0 0 0 255 255 255 255 255 255 255 255 255 255 255 0 +0 0 0 0 0 255 255 255 255 255 255 255 255 255 255 255 0 +0 0 0 0 0 0 0 255 255 255 255 255 255 255 0 0 0 +0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 diff --git a/skimage/morphology/data/two-hole_fiji.txt b/skimage/morphology/data/two-hole_fiji.txt new file mode 100644 index 00000000..6c16d3f7 --- /dev/null +++ b/skimage/morphology/data/two-hole_fiji.txt @@ -0,0 +1,16 @@ +0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 255 255 0 255 255 0 0 0 0 0 +0 0 0 0 0 255 255 0 0 255 0 0 255 0 0 0 0 +0 0 0 0 255 0 0 0 0 255 0 0 0 255 0 0 0 +0 0 0 255 0 0 0 0 0 255 0 0 0 0 255 0 0 +0 0 0 255 0 0 0 0 0 255 0 0 0 0 255 0 0 +0 0 0 0 255 0 0 0 0 255 0 0 0 0 255 0 0 +0 0 0 0 255 0 0 0 0 255 0 0 0 0 255 0 0 +0 0 0 0 0 255 0 0 0 255 0 0 0 0 255 0 0 +0 0 0 0 0 0 255 255 0 255 0 0 0 255 0 0 0 +0 0 0 0 0 0 0 0 255 0 255 255 255 0 0 0 0 +0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 diff --git a/skimage/morphology/skel.py b/skimage/morphology/skel.py new file mode 100644 index 00000000..e53c7ee6 --- /dev/null +++ b/skimage/morphology/skel.py @@ -0,0 +1,678 @@ +from __future__ import division, print_function, absolute_import + +import numpy as np + + +def _prepare_image(img_in): + """Convert to a binary image, pad the it w/ zeros, and ensure it's 3D. + """ + if img_in.ndim < 2 or img_in.ndim > 3: + raise ValueError('expect 2D, got ndim = %s' % img_in.ndim) + + img = img_in.copy() + + if img.ndim == 2: + img = img.reshape((1,) + img.shape) + + # normalize to binary + img[img != 0] = 1 + + # pad w/ zeros to simplify dealing w/ neighborhood of a pixel + img_o = np.zeros(tuple(s + 2 for s in img.shape), + dtype=np.uint8) + img_o[1:-1, 1:-1, 1:-1] = img.astype(np.uint8) + + return img_o + + +def _postprocess_image(img_o): + """Clip the image (padding is an implementation detail), convert to b/w. + """ + img_oo = img_o[1:-1, 1:-1, 1:-1] + img_oo = img_oo.squeeze() + img_oo *= 255 + return img_oo + + +def get_neighborhood(img, p, r, c): + """Get the neighborhood of a pixel. + + Assume zero boundary conditions. Image is already padded, so no + out-of-bounds checking. + """ + neighborhood = np.zeros(27, dtype=np.uint8) + + neighborhood[0] = img[p-1, r-1, c-1] + neighborhood[1] = img[p-1, r, c-1] + neighborhood[2] = img[p-1, r+1, c-1] + + neighborhood[ 3] = img[p-1, r-1, c] + neighborhood[ 4] = img[p-1, r, c] + neighborhood[ 5] = img[p-1, r+1, c] + + neighborhood[ 6] = img[p-1, r-1, c+1] + neighborhood[ 7] = img[p-1, r, c+1] + neighborhood[ 8] = img[p-1, r+1, c+1] + + neighborhood[ 9] = img[p, r-1, c-1] + neighborhood[10] = img[p, r, c-1] + neighborhood[11] = img[p, r+1, c-1] + + neighborhood[12] = img[p, r-1, c] + neighborhood[13] = img[p, r, c] + neighborhood[14] = img[p, r+1, c] + + neighborhood[15] = img[p, r-1, c+1] + neighborhood[16] = img[p, r, c+1] + neighborhood[17] = img[p, r+1, c+1] + + neighborhood[18] = img[p+1, r-1, c-1] + neighborhood[19] = img[p+1, r, c-1] + neighborhood[20] = img[p+1, r+1, c-1] + + neighborhood[21] = img[p+1, r-1, c] + neighborhood[22] = img[p+1, r, c] + neighborhood[23] = img[p+1, r+1, c] + + neighborhood[24] = img[p+1, r-1, c+1] + neighborhood[25] = img[p+1, r, c+1] + neighborhood[26] = img[p+1, r+1, c+1] + + return neighborhood + + +###### look-up tables +def fill_numpoints_LUT(n=256): + p = int(np.log2(n) + 1) + return np.sum(np.arange(n)[:, None] & (1 << np.arange(p)) != 0, axis=1) + +NUMPOINTS_LUT = fill_numpoints_LUT() + + +def fill_Euler_LUT(): + LUT = np.zeros(256, dtype=int) + + LUT[1] = 1 + LUT[3] = -1 + LUT[5] = -1 + LUT[7] = 1 + LUT[9] = -3 + LUT[11] = -1 + LUT[13] = -1 + LUT[15] = 1 + LUT[17] = -1 + LUT[19] = 1 + LUT[21] = 1 + LUT[23] = -1 + LUT[25] = 3 + LUT[27] = 1 + LUT[29] = 1 + LUT[31] = -1 + LUT[33] = -3 + LUT[35] = -1 + LUT[37] = 3 + LUT[39] = 1 + LUT[41] = 1 + LUT[43] = -1 + LUT[45] = 3 + LUT[47] = 1 + LUT[49] = -1 + LUT[51] = 1 + + LUT[53] = 1 + LUT[55] = -1 + LUT[57] = 3 + LUT[59] = 1 + LUT[61] = 1 + LUT[63] = -1 + LUT[65] = -3 + LUT[67] = 3 + LUT[69] = -1 + LUT[71] = 1 + LUT[73] = 1 + LUT[75] = 3 + LUT[77] = -1 + LUT[79] = 1 + LUT[81] = -1 + LUT[83] = 1 + LUT[85] = 1 + LUT[87] = -1 + LUT[89] = 3 + LUT[91] = 1 + LUT[93] = 1 + LUT[95] = -1 + LUT[97] = 1 + LUT[99] = 3 + LUT[101] = 3 + LUT[103] = 1 + + LUT[105] = 5 + LUT[107] = 3 + LUT[109] = 3 + LUT[111] = 1 + LUT[113] = -1 + LUT[115] = 1 + LUT[117] = 1 + LUT[119] = -1 + LUT[121] = 3 + LUT[123] = 1 + LUT[125] = 1 + LUT[127] = -1 + LUT[129] = -7 + LUT[131] = -1 + LUT[133] = -1 + LUT[135] = 1 + LUT[137] = -3 + LUT[139] = -1 + LUT[141] = -1 + LUT[143] = 1 + LUT[145] = -1 + LUT[147] = 1 + LUT[149] = 1 + LUT[151] = -1 + LUT[153] = 3 + LUT[155] = 1 + + LUT[157] = 1 + LUT[159] = -1 + LUT[161] = -3 + LUT[163] = -1 + LUT[165] = 3 + LUT[167] = 1 + LUT[169] = 1 + LUT[171] = -1 + LUT[173] = 3 + LUT[175] = 1 + LUT[177] = -1 + LUT[179] = 1 + LUT[181] = 1 + LUT[183] = -1 + LUT[185] = 3 + LUT[187] = 1 + LUT[189] = 1 + LUT[191] = -1 + LUT[193] = -3 + LUT[195] = 3 + LUT[197] = -1 + LUT[199] = 1 + LUT[201] = 1 + LUT[203] = 3 + LUT[205] = -1 + LUT[207] = 1 + + LUT[209] = -1 + LUT[211] = 1 + LUT[213] = 1 + LUT[215] = -1 + LUT[217] = 3 + LUT[219] = 1 + LUT[221] = 1 + LUT[223] = -1 + LUT[225] = 1 + LUT[227] = 3 + LUT[229] = 3 + LUT[231] = 1 + LUT[233] = 5 + LUT[235] = 3 + LUT[237] = 3 + LUT[239] = 1 + LUT[241] = -1 + LUT[243] = 1 + LUT[245] = 1 + LUT[247] = -1 + LUT[249] = 3 + LUT[251] = 1 + LUT[253] = 1 + LUT[255] = -1 + return LUT + +LUT = fill_Euler_LUT() + + +### Octants (indexOctantXXX functions) +OCTANTS = tuple(range(8)) +NEB, NWB, SEB, SWB, NEU, NWU, SEU, SWU = OCTANTS + +neib_idx = np.empty((8, 7), dtype=int) +neib_idx[NEB, ...] = [2, 1, 11, 10, 5, 4, 14] +neib_idx[NWB, ...] = [0, 9, 3, 12, 1, 10, 4] +neib_idx[SEB, ...] = [8, 7, 17, 16, 5, 4, 14] +neib_idx[SWB, ...] = [6, 15, 7, 16, 3, 12, 4] +neib_idx[NEU, ...] = [20, 23, 19, 22, 11, 14, 10] +neib_idx[NWU, ...] = [18, 21, 9, 12, 19, 22, 10] +neib_idx[SEU, ...] = [26, 23, 17, 14, 25, 22, 16] +neib_idx[SWU, ...] = [24, 25, 15, 16, 21, 22, 12] + +def index_octants(octant, neighbors): + n = 1 + for j, idx in enumerate(neib_idx[octant]): + if neighbors[idx] == 1: + n |= 2**(7 - j) + return n + + +def is_surfacepoint(neighbors, points_LUT): + for octant in OCTANTS: + n = index_octants(octabt, neighbors) + if n not in (240, 165, 170) and points_LUT[n] > 2: + return False + return True + + +def is_Euler_invariant(neighbors): + """Check if a point is Euler invariant. + + Calculate Euler characteristc for each octant and sum up. + + Parameters + ---------- + neighbors : ndarray, shape (27,) + neighbors of a point + + Returns + ------- + bool + + """ + euler_char = 0 + for octant in OCTANTS: + n = index_octants(octant, neighbors) + euler_char += LUT[n] + return euler_char == 0 + + +def is_simple_point(neighbors): + """Check is a point is a Simple Point. + + This method is named 'N(v)_labeling' in [Lee94]. + Outputs the number of connected objects in a neighborhood of a point + after this point would have been removed. + + Parameters + ---------- + neighbors : ndarray, shape(27,) + neighbors of the point + + Returns + ------- + bool + Whether the point is simple or not. + + """ + # copy neighbors for labeling + # ignore center pixel (i=13) when counting (see [Lee94]) + cube = np.r_[neighbors[:13], neighbors[14:]] + + # set initial label + label = 2 + + # for all point in the neighborhood + for i in range(26): + if cube[i] == 1: + # voxel has not been labeled yet + # start recursion with any octant that contains the point i + if i in (0, 1, 3, 4, 9, 10, 12): + octree_labeling(1, label, cube) + elif i in (2, 5, 11, 13): + octree_labeling(2, label, cube) + elif i in (6, 7, 14, 15): + octree_labeling(3, label, cube) + elif i in (8, 16): + octree_labeling(4, label, cube) + elif i in (17, 18, 20, 21): + octree_labeling(5, label, cube) + elif i in (19, 22): + octree_labeling(6, label, cube) + elif i in (23, 24): + octree_labeling(7, label, cube) + elif i == 25: + octree_labeling(8, label, cube) + else: + raise ValueError("Never be here. i = %s" % i) + label += 1 + if label - 2 >= 2: + return False + return True + + +def octree_labeling(octant, label, cube): + """This is a recursive method that calculates the number of connected + components in the 3D neighborhood after the center pixel would + have been removed. + + Parameters + ---------- + octant : int + octant index + label : int + the current label of the center point + cube : ndarray, shape(26,) + local neighborhood of the point + + """ + # check if there are points in the octant with value 1 + if octant == 1: + # set points in this octant to current label + # and recursive labeling of adjacent octants + if cube[0] == 1: + cube[0] = label + if cube[1] == 1: + cube[1] = label + octree_labeling(2, label, cube) + if cube[3] == 1: + cube[3] = label + octree_labeling(3, label, cube) + if cube[4] == 1: + cube[4] = label + octree_labeling(2, label, cube) + octree_labeling(3, label, cube) + octree_labeling(4, label, cube) + if cube[9] == 1: + cube[9] = label + octree_labeling(5, label, cube) + if cube[10] == 1: + cube[10] = label + octree_labeling(2, label, cube) + octree_labeling(5, label, cube) + octree_labeling(6, label, cube) + if cube[12] == 1: + cube[12] = label + octree_labeling(3, label, cube) + octree_labeling(5, label, cube) + octree_labeling(7, label, cube) + + if octant == 2: + if cube[1] == 1: + cube[1] = label + octree_labeling(1, label, cube) + if cube[4] == 1: + cube[4] = label + octree_labeling(1, label, cube) + octree_labeling(3, label, cube) + octree_labeling(4, label, cube) + if cube[10] == 1: + cube[10] = label + octree_labeling(1, label, cube) + octree_labeling(5, label, cube) + octree_labeling(6, label, cube) + if cube[2] == 1: + cube[2] = label + if cube[5] == 1: + cube[5] = label + octree_labeling(4, label, cube) + if cube[11] == 1: + cube[11] = label + octree_labeling(6, label, cube) + if cube[13] == 1: + cube[13] = label + octree_labeling(4, label, cube) + octree_labeling(6, label, cube) + octree_labeling(8, label, cube) + + if octant ==3: + if cube[3] == 1: + cube[3] = label + octree_labeling(1, label, cube) + if cube[4] == 1: + cube[4] = label + octree_labeling(1, label, cube) + octree_labeling(2, label, cube) + octree_labeling(4, label, cube) + if cube[12] == 1: + cube[12] = label + octree_labeling(1, label, cube) + octree_labeling(5, label, cube) + octree_labeling(7, label, cube) + if cube[6] == 1: + cube[6] = label + if cube[7] == 1: + cube[7] = label + octree_labeling(4, label, cube) + if cube[14] == 1: + cube[14] = label + octree_labeling(7, label, cube) + if cube[15] == 1: + cube[15] = label + octree_labeling(4, label, cube) + octree_labeling(7, label, cube) + octree_labeling(8, label, cube) + + if octant == 4: + if cube[4] == 1: + cube[4] = label + octree_labeling(1, label, cube) + octree_labeling(2, label, cube) + octree_labeling(3, label, cube) + if cube[5] == 1: + cube[5] = label + octree_labeling(2, label, cube) + if cube[13] == 1: + cube[13] = label + octree_labeling(2, label, cube) + octree_labeling(6, label, cube) + octree_labeling(8, label, cube) + if cube[7] == 1: + cube[7] = label + octree_labeling(3, label, cube) + if cube[15] == 1: + cube[15] = label + octree_labeling(3, label, cube) + octree_labeling(7, label, cube) + octree_labeling(8, label, cube) + if cube[8] == 1: + cube[8] = label + if cube[16] == 1: + cube[16] = label + octree_labeling(8, label, cube) + + if octant == 5: + if cube[9] == 1: + cube[9] = label + octree_labeling(1, label, cube) + if cube[10] == 1: + cube[10] = label + octree_labeling(1, label, cube) + octree_labeling(2, label, cube) + octree_labeling(6, label, cube) + if cube[12] == 1: + cube[12] = label + octree_labeling(1, label, cube) + octree_labeling(3, label, cube) + octree_labeling(7, label, cube) + if cube[17] == 1: + cube[17] = label + if cube[18] == 1: + cube[18] = label + octree_labeling(6, label, cube) + if cube[20] == 1: + cube[20] = label + octree_labeling(7, label, cube) + if cube[21] == 1: + cube[21] = label + octree_labeling(6, label, cube) + octree_labeling(7, label, cube) + octree_labeling(8, label, cube) + + if octant == 6: + if cube[10] == 1: + cube[10] = label + octree_labeling(1, label, cube) + octree_labeling(2, label, cube) + octree_labeling(5, label, cube) + if cube[11] == 1: + cube[11] = label + octree_labeling(2, label, cube) + if cube[13] == 1: + cube[13] = label + octree_labeling(2, label, cube) + octree_labeling(4, label, cube) + octree_labeling(8, label, cube) + if cube[18] == 1: + cube[18] = label + octree_labeling(5, label, cube) + if cube[21] == 1: + cube[21] = label + octree_labeling(5, label, cube) + octree_labeling(7, label, cube) + octree_labeling(8, label, cube) + if cube[19] == 1: + cube[19] = label + if cube[22] == 1: + cube[22] = label + octree_labeling(8, label, cube) + + if octant == 7: + if cube[12] == 1: + cube[12] = label + octree_labeling(1, label, cube) + octree_labeling(3, label, cube) + octree_labeling(5, label, cube) + if cube[14] == 1: + cube[14] = label + octree_labeling(3, label, cube) + if cube[15] == 1: + cube[15] = label + octree_labeling(3, label, cube) + octree_labeling(4, label, cube) + octree_labeling(8, label, cube) + if cube[20] == 1: + cube[20] = label + octree_labeling(5, label, cube) + if cube[21] == 1: + cube[21] = label + octree_labeling(5, label, cube) + octree_labeling(6, label, cube) + octree_labeling(8, label, cube) + if cube[23] == 1: + cube[23] = label + if cube[24] == 1: + cube[24] = label + octree_labeling(8, label, cube) + + if octant == 8: + if cube[13] == 1: + cube[13] = label + octree_labeling(2, label, cube) + octree_labeling(4, label, cube) + octree_labeling(6, label, cube) + if cube[15] == 1: + cube[15] = label + octree_labeling(3, label, cube) + octree_labeling(4, label, cube) + octree_labeling(7, label, cube) + if cube[16] == 1: + cube[16] = label + octree_labeling(4, label, cube) + if cube[21] == 1: + cube[21] = label + octree_labeling(5, label, cube) + octree_labeling(6, label, cube) + octree_labeling(7, label, cube) + if cube[22] == 1: + cube[22] = label + octree_labeling(6, label, cube) + if cube[24] == 1: + cube[24] = label + octree_labeling(7, label, cube) + if cube[25] == 1: + cube[25] = label + + +def _loop_through(img, curr_border): + """Inner loop of compute_thin_image. + + return simple_border_points as a list to be rechecked sequentially. + """ + # loop through the image + # NB: each loop is from 1 to size-1: img is padded from all sides + simple_border_points = [] + + ### XXX: 2D images + ### if the original is 2D, img.shape[0] == 3, the algorithm removes too much + ### because all points are considered 'boundary' in the 3rd direction. + ### Hence just bail out + if img.shape[0] == 3 and curr_border in (5, 6): + print("skipping curr_border = ", curr_border) + return [] + + for p in range(1, img.shape[0] - 1): + for r in range(1, img.shape[1] - 1): + for c in range(1, img.shape[2] - 1): + + # check if pixel is foreground + if img[p, r, c] != 1: + continue + + is_border_pt = (curr_border == 1 and img[p, r, c-1] <= 0 or #N + curr_border == 2 and img[p, r, c+1] <= 0 or #S + curr_border == 3 and img[p, r+1, c] <= 0 or #E + curr_border == 4 and img[p, r-1, c] <= 0 or #W + curr_border == 5 and img[p+1, r, c] <= 0 or #U + curr_border == 6 and img[p-1, r, c] <= 0) #B + if not is_border_pt: + # current point is not deletable + continue + + neighborhood = get_neighborhood(img, p, r, c) + + # check if (p, r, c) is an endpoint. An endpoint has exactly + # one neighbor in the 26-neighborhood. + # The center pixel is counted, thus r.h.s. is 2 + if neighborhood.sum() == 2: + continue + + # check if point is Euler invariant (condition 1 in [Lee94]) + # if it is not, it's not deletable + if not is_Euler_invariant(neighborhood): + continue + + # check if point is simple (i.e., deletion does not + # change connectivity in the 3x3x3 neighborhood) + # this are conditions 2 and 3 in [Lee94] + if not is_simple_point(neighborhood): + continue + + # ok, add (p, r, c) to the list of simple border points + simple_border_points.append((p, r, c)) + return simple_border_points + + +def compute_thin_image(img_in): + + ### prepare + img = _prepare_image(img_in) + + ### compute + unchanged_borders = 0 + + # loop through the image several times until there is no change for all + # the six border types + while unchanged_borders < 6: + unchanged_borders = 0 + for curr_border in (4, 3, 2, 1, 5, 6): + + simple_border_points = _loop_through(img, curr_border) + print(curr_border, " : ", simple_border_points, '\n') + + # sequential re-checking to preserve connectivity when deleting + # in a parallel way + no_change = True + for pt in simple_border_points: + p, r, c = pt + neighb = get_neighborhood(img, p, r, c) + if is_simple_point(neighb): + img[p, r, c] = 0 + no_change = False + else: + print(" *** ", pt, is_simple_point(neighb)) + + if no_change: + unchanged_borders += 1 + simple_border_points = [] + + img = _postprocess_image(img) + return img + + +if __name__ == "__main__": + pass diff --git a/skimage/morphology/test_simple_2d.py b/skimage/morphology/test_simple_2d.py new file mode 100644 index 00000000..c4a5cd2b --- /dev/null +++ b/skimage/morphology/test_simple_2d.py @@ -0,0 +1,84 @@ +from __future__ import division, print_function, absolute_import + +import numpy as np +from numpy.testing import assert_equal + +from skimage import io + +from skel import compute_thin_image + +import matplotlib.pyplot as plt +import matplotlib.ticker as ticker + +# nose test generators: + +# 2D images +def test_simple_2d_images(): + for fname in ("strip", "loop", "cross", "two-hole"): + yield check_skel, fname + +# trivial 3D images +def test_simple_3d(): + for fname in ['3/stack', '4/stack']: + yield check_skel_3d, fname + + +def check_skel(fname, viz=False): + # compute the thin image and compare the result to that of ImageJ + img = np.loadtxt('data/' + fname + '.txt', dtype=np.uint8) + + if viz: + ax = _viz(img, **dict(marker='s', color='b', s=99, alpha=0.2)) + + # compute + img1_2d = compute_thin_image(img) + + if viz: + ax = _viz(img1_2d, ax, **dict(marker='o', color='r', + s=80, alpha=0.7, label='us')) + + # compare to FIJI + img_f = np.loadtxt('data/' + fname + '_fiji.txt', dtype=np.uint8) + + if not viz: + # actually compare images + assert_equal(img1_2d, img_f) + else: + ax = _viz(img_f, ax, **dict(marker='o', color='g', s=45, label='fiji')) + + ax.legend() + ax.grid(True) + + def yformatter(val, pos): + return int(img.shape[1] - val + 1) + def xformatter(val, pos): + return int(val + 1) + ax.xaxis.set_major_formatter(ticker.FuncFormatter(xformatter)) + ax.yaxis.set_major_formatter(ticker.FuncFormatter(yformatter)) + + plt.show() + + +def _viz(img, ax=None, **kwds): + if ax is None: + import matplotlib.pyplot as plt + fix, ax = plt.subplots() + + x, y = np.nonzero(img) + ax.scatter(y, img.shape[1] - x, **kwds) + return ax + + +def check_skel_3d(fname): + img = io.imread('data/' + fname + '.tif') + img_f = io.imread('data/' + fname + '_fiji.tif') + + img_s = compute_thin_image(img) + assert_equal(img_s, img_f) + + +if __name__ == "__main__": + import sys + if len(sys.argv) < 2: + sys.exit("Expect an image name from the data/ directory.") + check_skel(sys.argv[1], True) diff --git a/skimage/morphology/trivial.py b/skimage/morphology/trivial.py new file mode 100644 index 00000000..46656d93 --- /dev/null +++ b/skimage/morphology/trivial.py @@ -0,0 +1,73 @@ +import numpy as np +import matplotlib.pyplot as plt + +from skel import prepare_image, compute_thin_image + + +def get_rhombus(n=64, L=22, width=7): + img = np.zeros((n, n), dtype=np.int8) + + x = np.arange(L, dtype=int) + y = L - x + + for w in range(width): + img[x + n//2, y + n//2 + w] = 1 + img[-x + n//2, y + n//2 + w] = 1 + img[x + n//2, -y + n//2 + w] = 1 + img[-x + n//2, -y + n//2 + w] = 1 + + return img + +def get_strip(): + img = np.array([[0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0], + [0, 0, 0, 0, 1, 1, 1, 0, 0, 0, 0, 0], + [0, 0, 0, 0, 1, 1, 1, 0, 0, 0, 0, 0], + [0, 0, 0, 0, 1, 1, 1, 0, 0, 0, 0, 0], + [0, 0, 0, 0, 1, 1, 1, 0, 0, 0, 0, 0], + [0, 0, 0, 0, 1, 1, 1, 0, 0, 0, 0, 0], + [0, 0, 0, 0, 1, 1, 1, 0, 0, 0, 0, 0], + [0, 0, 0, 0, 1, 1, 1, 0, 0, 0, 0, 0], + [0, 0, 0, 0, 1, 1, 1, 0, 0, 0, 0, 0], + [0, 0, 0, 0, 1, 1, 1, 0, 0, 0, 0, 0], + [0, 0, 0, 0, 1, 1, 1, 0, 0, 0, 0, 0], + [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]], dtype=np.uint8) + return img + + +def get_square(n=64, L=16, width=4): + img = np.zeros((n, n), dtype=np.int8) + + x = np.arange(-L//2 + 1, L//2) + for w in range(width): + img[n//2 + x, n//2 + L//2 + w] = 1 + img[n//2 + x, n//2 - L//2 - w] = 1 + + img[n//2 + L//2 - w, n//2 + x] = 1 + img[n//2 + x, n//2 - L//2 - w] = 1 + + return img + + +def get_loop(): + img = np.loadtxt('Untitled.txt', dtype=np.uint8) + return img + + +if __name__ == "__main__": + +## img = get_rhombus() +## img = get_strip() + img = get_loop() + + x, y = np.nonzero(img) + plt.scatter(x, y, marker='s', color='b', s=40, alpha=0.3) + + # skeletonize + img1 = prepare_image(img) + img1 = compute_thin_image(img1) + + img1_2d = img1[1, 1:, 1:] + x, y = np.nonzero(img1_2d) + + plt.scatter(x, y, marker='o', color='r') + plt.show()