From e5dab80e9c5bf1d4f197ffc820ae07ab7a33e910 Mon Sep 17 00:00:00 2001 From: Chintak Sheth Date: Tue, 23 Apr 2013 19:47:10 +0530 Subject: [PATCH] Clean, and implemented using horse.jpeg --- doc/examples/applications/plot_morphology.py | 272 +++++++++++++++++++ skimage/data/horse.jpeg | Bin 0 -> 23528 bytes 2 files changed, 272 insertions(+) create mode 100644 doc/examples/applications/plot_morphology.py create mode 100644 skimage/data/horse.jpeg diff --git a/doc/examples/applications/plot_morphology.py b/doc/examples/applications/plot_morphology.py new file mode 100644 index 00000000..206375ef --- /dev/null +++ b/doc/examples/applications/plot_morphology.py @@ -0,0 +1,272 @@ +""" +======================= +Morphological Filtering +======================= + +Morphological image processing is a collection of non-linear operations related +to the shape or morphology of features in an image, such as boundaries, +skeletons, etc. In any given technique, we probe an image with a small shape or +template called a structuring element, which defines the region of interest or +neighborhood around a pixel. + +In this document we outline the following basic morphological operations: + +1. Erosion +2. Dilation +3. Opening +4. Closing +5. White Tophat +6. Black Tophat +7. Skeletonize +8. Convex Hull + + +To get started, let's load an image using ``io.imread``. Note that morphology +functions only work on gray-scale or binary images, so we set ``as_grey=True``. +""" + +import matplotlib.pyplot as plt +from skimage.data import data_dir +from skimage.util import img_as_ubyte +from skimage import io + +plt.gray() +phantom = img_as_ubyte(io.imread(data_dir+'/phantom.png', as_grey=True)) +plt.imshow(phantom) + +""" +.. image:: PLOT2RST.current_figure + +Let's also define a convenience function for plotting comparisons: +""" + +def plot_comparison(original, filtered, filter_name): + + fig, (ax1, ax2) = plt.subplots(ncols=2, figsize=(8, 4)) + ax1.imshow(original) + ax1.set_title('original') + ax1.axis('off') + ax2.imshow(filtered) + ax2.set_title(filter_name) + ax2.axis('off') + +""" +Erosion +======= + +Morphological ``erosion`` sets a pixel at (i, j) to the *minimum over all +pixels in the neighborhood centered at (i, j)*. The structuring element, +``selem``, passed to ``erosion`` is a boolean array that describes this +neighborhood. Below, we use ``disk`` to create a circular structuring element, +which we use for most of the following examples. +""" + +from skimage.morphology import erosion, dilation, opening, closing, white_tophat +from skimage.morphology import black_tophat, skeletonize, convex_hull_image +from skimage.morphology import disk + +selem = disk(6) +eroded = erosion(phantom, selem) +plot_comparison(phantom, eroded, 'erosion') + +""" +.. image:: PLOT2RST.current_figure + +Notice how the white boundary of the image disappears or gets eroded as we +increase the size of the disk. Also notice the increase in size of the two +black ellipses in the center and the disappearance of the 3 light grey +patches in the lower part of the image. + + +Dilation +======== + +Morphological ``dilation`` sets a pixel at (i, j) to the *maximum over all +pixels in the neighborhood centered at (i, j)*. Dilation enlarges bright +regions and shrinks dark regions. +""" + +dilated = dilation(phantom, selem) +plot_comparison(phantom, dilated, 'dilation') + +""" +.. image:: PLOT2RST.current_figure + +Notice how the white boundary of the image thickens, or gets dilated, as we +increase the size of the disk. Also notice the decrease in size of the two +black ellipses in the centre, and the thickening of the light grey circle in +the center and the 3 patches in the lower part of the image. + + +Opening +======= + +Morphological ``opening`` on an image is defined as an *erosion followed by a +dilation*. Opening can remove small bright spots (i.e. "salt") and connect +small dark cracks. +""" + +opened = opening(phantom, selem) +plot_comparison(phantom, opened, 'opening') + +""" +.. image:: PLOT2RST.current_figure + +Since ``opening`` an image starts with an erosion operation, light regions that +are *smaller* than the structuring element are removed. The dilation operation +that follows ensures that light regions that are *larger* than the structuring +element retain their original size. Notice how the light and dark shapes in the +center their original thickness but the 3 lighter patches in the bottom get +completely eroded. The size dependence is highlighted by the outer white ring: +The parts of the ring thinner than the structuring element were completely +erased, while the thicker region at the top retains its original thickness. + + +Closing +======= + +Morphological ``closing`` on an image is defined as a *dilation followed by an +erosion*. Closing can remove small dark spots (i.e. "pepper") and connect +small bright cracks. + +To illustrate this more clearly, let's add a small crack to the white border: +""" + +phantom = img_as_ubyte(io.imread(data_dir+'/phantom.png', as_grey=True)) +phantom[10:30, 200:210] = 0 + +closed = closing(phantom, selem) +plot_comparison(phantom, closed, 'closing') + +""" +.. image:: PLOT2RST.current_figure + +Since ``closing`` an image starts with an dilation operation, dark regions +that are *smaller* than the structuring element are removed. The dilation +operation that follows ensures that dark regions that are *larger* than the +structuring element retain their original size. Notice how the white ellipses +at the bottom get connected because of dilation, but other dark region retain +their original sizes. Also notice how the crack we added is mostly removed. + + +White tophat +============ + +The ``white_tophat`` of an image is defined as the *image minus its +morphological opening*. This operation returns the bright spots of the image +that are smaller than the structuring element. + +To make things interesting, we'll add bright and dark spots to the image: +""" + +phantom = img_as_ubyte(io.imread(data_dir+'/phantom.png', as_grey=True)) +phantom[340:350, 200:210] = 255 +phantom[100:110, 200:210] = 0 + +w_tophat = white_tophat(phantom, selem) +plot_comparison(phantom, w_tophat, 'white tophat') + +""" +.. image:: PLOT2RST.current_figure + +As you can see, the 10-pixel wide white square is highlighted since it is +smaller than the structuring element. Also, the thin, white edges around most +of the ellipse are retained because they're smaller than the structuring +element, but the thicker region at the top disappears. + + +Black tophat +============ + +The ``black_tophat`` of an image is defined as its morphological **closing +minus the original image**. This operation returns the *dark spots of the +image that are smaller than the structuring element*. +""" + +b_tophat = black_tophat(phantom, selem) +plot_comparison(phantom, b_tophat, 'black tophat') + +""" +.. image:: PLOT2RST.current_figure + +As you can see, the 10-pixel wide black square is highlighted since it is +smaller than the structuring element. + + +Duality +------- + +As you should have noticed, many of these operations are simply the reverse +of another operation. This duality can be summarized as follows: + +1. Erosion <-> Dilation +2. Opening <-> Closing +3. White tophat <-> Black tophat + + +Skeletonize +=========== + +Thinning is used to reduce each connected component in a binary image to a +*single-pixel wide skeleton*. It is important to note that this is performed +on binary images only. + +""" + +from skimage import img_as_bool +horse = ~img_as_bool(io.imread(data_dir+'/horse.jpeg', as_grey=True)) + +sk = skeletonize(horse) +plot_comparison(horse, sk, 'skeletonize') + +""" +.. image:: PLOT2RST.current_figure + +As the name suggests, this technique is used to thin the image to 1-pixel wide +skeleton by applying thinning successively. + + +Convex hull +=========== + +The ``convex_hull_image`` is the *set of pixels included in the smallest +convex polygon that surround all white pixels in the input image*. Again note +that this is also performed on binary images. + +""" + +hull1 = convex_hull_image(horse) +plot_comparison(horse, hull1, 'convex hull') + +""" +.. image:: PLOT2RST.current_figure + +As the figure illustrates, ``convex_hull_image`` gives the smallest polygon +which covers the white or True completely in the image. + +If we add a small grain to the image, we can see how the convex hull adapts to +enclose that grain: +""" + +import numpy as np + +horse2 = np.copy(horse) +horse2[45:50, 75:80] = 1 + +hull2 = convex_hull_image(horse2) +plot_comparison(horse2, hull2, 'convex hull') + +""" +.. image:: PLOT2RST.current_figure + + +Additional Resources +==================== + +1. `MathWorks tutorial on morphological processing `_ +2. Auckland university's tutorial on Morphological Image Processing http://www.cs.auckland.ac.nz/courses/compsci773s1c/lectures/ImageProcessing-html/topic4.htm +3. http://en.wikipedia.org/wiki/Mathematical_morphology + +""" + +plt.show() \ No newline at end of file diff --git a/skimage/data/horse.jpeg b/skimage/data/horse.jpeg new file mode 100644 index 0000000000000000000000000000000000000000..f37e1a312a9408e0f0222547ac56103fea1240bf GIT binary patch literal 23528 zcmd43d0Z3OwlVOrnN(zwib^O=rTDhp=bS$0-0pktd*1ht4}X#nD#aewTFMYW?)KYzet(>iIpgnF5HDmZ{p0!9-x8+#365XK z@PGC7^(OyUH~9MI=~{x;EaDHu3ra*S!VE2gmrnSh~#go8{kHTCKBQzrn`w`%O;HKWyHzd(Yl|uKN$Tc^&aS z>f`I@9}*gNBAgj#&n_C|G;(np3mpeRuIjVV$bU+Fgjw;o~pVbeC9) zmP)@)?YEiz?Oce8OxcBG^17EdoTq@oxy3GosZ9|Nnhx>HH#5^+Yp0xrQQl zN)rqCE9zl&Y-gZUU3}uKDBwcOG@*g3++UkPmfGrZ&GNIP#*ZIKRUaaI2j0i#6QW$E z2{UcL&1=x!CJJ(ZjG9_tb9pjHDSqP)AN7&(P?DJ2Z{6fMyL;$7V~!(3Vn;c^Q)Ikn z6t~{`_{q+&_NP(fg=oobtw`^?{st)tRU7VB-#~bYmI2Snqox1vHb^NBjgzhr3F;wL z|3N!6PpVW0&#uzA(Y&mu%126SvZ2JXmS;7^D1YIcVL3+}x-~4kK5m_ICR;hgo)XP1 zF2jtRrwPg@6zTb48YWI?zr=#?(wKVWgF{3XyG=l5VD3LIRZ1w7P7{(XsK%^Sy=MC_ z;#WjG!+8cP9Q@KeT`AskC16X z9My=meOKd$K)uF>iQvfQS)^wC+i61VD_+2gCkvP+vgL~OPqRLrA)TuVk*;kv=T5$Q z*GpB^9(rl<)>r960-jKlL{b`l%yTt>B4!x7H znd^!I?!F6%JZ`x6c}dz_ek1|@o{soIV|fxD9}GENY|&3{PH>kxCqUi;HB=&o&S8<_ zSevCSlCXp`wi+ULJD+`BZvxuD%0IhSfc=1;{*u=}IFL9^P*qfXj#Eq~ULa2z z{XcSSO{jTrJwk7wtD>{2dM@HZNyM#DNbWG4qK!YL7oR8-divZLc$4^CDC*m+9P8qu+h3M4*AD*!P(I*kzd<;M5a3 z{C=#&Y|WRB{XxeHR_3SnjMI-7&eo3(>-v4De;e_C_CQI6klhyOw+Ul8z>OcE%h)UM z7oMh9hZ-{$8pLJSb7DQ`r9K-K`!uz&>*sX7M8mg~KXSP9&c?CLA2*YhKI~yrkYB_1 zU928>&%5jEru&eHyTZsG0csXLehoxHcWHGaPp@jvAOY356`*I3(7@{}s6ClJUe>^w zy-Dfea0T}qUAGLLi!m5MC!R1oghH(M^=LwTtro4suKCzrYccH%&mT;|{r>M}cfR^{(?Cz6&_Y$cQe9xoEbe>g*OQ{_ zb|T$qnsC;tT744_3~Uj%-&FRiB6FI+pmEr16t`eM?gCGANNkt=1chohL^3C(FR{n?0sFYsJ z@wd^;xf76THQszJ*LO;G^3LVEeXXwt_WQfBqz7Mm_CJHV4)O;YH0mGE`nL9Xc_w*k zl_C!EGqRA?DUx|d>HKh)s{mI;;Vc~#mJE3E+pY*ZSgyIulZ~jL;o3VNq&MkMZ*M=KCsJw3rDs;zNQ9}N@DGG{%g z2)fJ&vFgdXuQ=Z^)kHO|>apE?&q8*!D8ini@^5SXTYS0Az~-mRkGr#WD=y*Ycs+Y% z53AjP{p+Ab1feso_OVyQ{CSyJY!=5?_9Q(=0r-`~;YR-dI(r%!Ph6+t? zSZ;vZVBtLvos!PbHq&h-8L}_3@us}=+m#zv?br-Vba6XsA;5H98>b1#ULHTVc$zSS zKJvw4EF5v%-D63QW|PaUbc@{T0b__KLYvLWOLsDEld*l>TgfU5-i&+_hCis)X;brHhJsXDP+KQKAMf@Hk? z7W3_6xBfd}W{wTd)yWG`cB@($dnJ;cPhx1o@k;U3Zxyc;O3KGnET z%$SF!xpuPbL_AGqQ-!fKnVTaT4XTBI4Yt+2SVidUU6 z`1I+z@E!KcN8EyLWQ|yV*;n8de#bx67zp4!-@~@1B3F+ z%?5iubPPN4x6D646i2W5`;e3r*0Ase9sO2=H$Go(HUdA;Ya-9NhOfn+oOJESbvNG~ z%uN}zRHRKB?F{@hy7AC5H?c(czB-C@+WfPskFYviGG|xA!%7{$6cclDLd8*$~`4z0sC-vxB zn&rpB0n1+-i*z@v47r|Rzxfm$7nSz@G<_-k>2Plw&uDK-Bx{s#GW+?7wRjCirva!XAyP=}M&@7ds@eLBF zA$fNH!Dr)Xfi|uCQQF~5#dUwr6`xZl?V*cs|8870P|gAqvnyE^D=Sc>;tN^U$Yj&d zOUw2%vQim)4eFA&%`&(2DnGo5tl}BKo@JB#e44O@(=e6Vol$dc_c-4uO~s!0dHk!7 z5oMGy6S)_*O%Xj!2%yC9N)qZ_k_W#*H;eo2+(1OR^gbZaq8TT5-AMaE%8L%WnFYNZ zv_MMmJkjWA+-!d2K>tW;#5Sk;r&4=}KOsA+xXISXdtd^4F5UsFSQU!N@z8u*y)0Bx zZ$17x@o}IP6YdDLJ32TR+0dbPkBE=(ZCaqQJlRp4-9%Al**%C^L;W;GUD+EN=}o#mcQD!~DVn_i z?*M@_bOO^Fl7$93Zhuj;0nL_OmO5QO5>L~_pBWRoSswl$!_MZhqQs+I0_)qZyk9V$ zUxzK@?l?W|(yCg@sc8b>48key{PtkjOTHSf&zRqFOdR2qA}-)SI+JuQ5Bm3T*Fk@X zj1e=Lhb)_w-;A44myJ47K43EmY%k^%PM!s9h=*3^lM{lb2}fyhhY!hbAr^HTlWAkl z98|_QJGd4-O;glF6=gI>Q@y%j;X4j$?s@+s1r6kljYiND11)4>il4{#5ud9^UPc@y zY2Uect{%PIPFV0O1qmTxdNWj$r1+rmQuAl-#{??y(f(g}hLPM@MqsXU!~X@TsxBBf@P-vz?bN~Ec*k#>4;W%UM?edS*n$|?|9GIj=#F$(7I7M!aKc^JQ3c}vBtjiJhx2v&25%W~0HVw{o2tit#YwS(l z7`Qf6Z+_OH?{$eN=gC33U&nxjKT}GqUyj$z#Rn)ewm)J&Z}1U8XYi$5AvcM^FkNuK zZ%_|y@evyNxK)snE9|bwW-Xb!VIk$-I_sX3S8kWy$R2bbS!aErCBzT*3D)VpJy2Zh zGF#T!f3+;5@8>t>^DaDHlN&z6-tu1sngE;oMQaH;Q)k~%c|((#G8<`i$q?DUM$O9s=no1FaWy{VN_dLGx_Ytk$o+1pHF<^s16LEaGgzB($8RJE{hq; zFH1>Zm`|URq<~q?{F-K00Jm%B$Fw(6)ZxSM>^AuVmPz-nDVIv=i2nm7N5+wwU!xw< zAe{@$lYk&oPcP*wf5(6*@!SLEhne)+MP$0l@&)LJQCMl{ZU`Z{9iN=@}wUak~8dWf$o5i zE$N05m7+U|G=01mnrBOw9)BD8OOjiXrxMcCf`P#&@b}i$4rJT$PTM+2qBP1WUQ(&) zqId#y_Y2`0rC+49INQ~d0?hzpEkFl8k@cC*jj3US7sVE8XwIewQoV#Khg&0aK3lu( zE6N%t5HkK57k8z@iq40zE|rtt)>rJ*a> zEj4X!?{oZK?vw+?e|mX1QsLghn_;$pz~9xzD7m67qJXuW92K7B#^|C_R`~hW%zI%E z#h&|$aO8;qd^TWIF)_@(D!McWS1pnda{>lcrU|Fr4=A#5HF`wIKtGeTv$KRpX2Rl> z({H(&U)X0Qx>D)r{}S@Z9Wda8lZM&j*dkbcn0JLSbCs7NaFVAS!*jAf^O&YMWU=lQxC zz>SdByOf;Dz}~aAJCOqSqSw^uJ`HN9sX18#YaPzebp0a3C z=N1EL3iorMBD=uj>KAKEqJ*~OQqdEWqQ~Ra11~EpuYcHVYPWgGP-n)^M4f+Qs!U7VASYJd6#JuyiTGlBLDVeH-=v;b=X2gvv_F|VKpq~DEoO5 z(S&<*w{L5kO3Hu#c)@y7G`UL<$3-H`#+`31ag_LsIUw3%?fh*GEDJPIswbREJ*hQI z8d%Wcs1kW~jTaDwsIL(iGnUKkJe$;6_KH}5%>E&G4Mf7|2y}jN!UD=-i91N1-acxV z@4L#J$|MIp1l|Y8hA{oODYg|lYSF)aaO?3>wlQkNO6XL;&2nv2!m}i*WC6wP0bFO7 zm3;n++#K%u3w1GYMNfT>hp;TmsNkWKt9t$L3i?^t#D zhRDj+wCoGl^q-R|ey5^D+*5#@fUvnKb(BKRI7JBeD+w*)0;^0Bnyj>|YFnU$L-M^! z_W>(8P2lj?)a8!O^+o#Yj55yNVC{Q8aYyRO6^(8i^%ldLm}3C>fv>fjaoS-bb_Osb z$_P-0RqliBF{#=l*JMBQT&O6hu|?TY+Ol(syejJ(aJ`)2%IlpIKq{na>Wjq^k##N`} zme>!9d$6xrL@`Wc*>&#W51&xs)90Zdn@+|?80#UV7|YhQG>3a=K_ma7T<1z*u$R}7 zvW$hj1;rDZ?<#t-4ij;Ez{2*>(G}!zF!Y(Os0fyBlyJ|y4+RCCD35)ZwGeJAzeKCL3Rf^0=;G8|V0 zHVd@!2aqt=b0L~`2Y8yhkoBcJs$V-PPTJxT`hdwle(WPhIxVH-FEu{we$WU*tqkFeH2h_laOX zGP)8T^9D>=r5c%PurTx&P_WI5l3X;olsdA;h+vNGir}Ae4Oe-voH}f^Z(V<_(oMAW z`H3Vq_&fh4y-mS%vT;?3h=wkoCXCqu7QDCD%8;Ug=IpfzPp;8q4(&J>kVscz9YPZL z=C(%pg1qXyae-`Erg;ft7Fss8F1#Wdp`W#NP^8B03Rsg#c`HIJaE*fiQ0g%sG+DwU z1aMWV)Dg{wTD|bpu$D&Nc{Y^bf`*GLA0PEym{5C?(=cu#W}IPp1VjdERTlMF7O0md z98Fr>6=nzR==ZUwlyH%)FuHb1!fBxKBSk~aH1r$OMo#pX-g=8T>heyOUjPy=T6a%e zU?jly_hQyjDmQNUdRlhb9uPt^efgAzp04p@7FP1_>GS-Jzv;|l*sklbF$ggADEf-6 z|5)^C!b%3ZQ5)d%eRdD8Z0W?+&s2l*gS^l22l<00 zR&9+%TJu)b&7SW!boi(HM4u2&RRkEwCJKL=Al{UM8KwgVoF-tpiVHPHb_oe^5ZrIt z!2e{tx@*k)bO26~@;I!GhhX{l*!F)aAxKde}~^zLhqRhi^ZhyF|S1UFZD&FycHUF(Fr=r7x4*A;gsX%By3 z&%)+Qd(VGz-zY#Qdsv6Eeaq%BDGfy#0^}0&>5v2Rg*)Nz^7w5u)DYQNW71pYEw+Q3 z@a5eME1otp$GscdB(`JhOl#}1N{;lb;qL6Qa{v64V>2>M*!8Hqs_M9^PdewgJ;@I2 z2^%QAb!9+EXv+Rg?)Fn$$2N%VZ0y`cQ?>ocF0u?X**lob4;*@dlMnip;hIG?Wk=Ta zrg|1i1NNc8w*pRtCya%snrB5C=DihpFPuq3*K=j#JVg~OyelwfIot}cxW`(Cs@@l{ z+rU=lWnDYc!4u6#58d~7gU{J+&W?k;D~vay$xu?&0k51uy_AzGS)w+(*8Gd!(d7F7 zOeZ1BGy%dhc!yfsUwYxHbi~C$y6-&SxrxD7eLqd`_4Hro^QS+NMw~INaM|V#D98d`J~al*H=y^M&34yjEI~7BR-v3b9?jkly72|9 zWiPBa=kRsN8(c7{J7^B(5*tQV3tB=h8NFJ2=-P22-tB+ZAUF8#_j~zw@2^{E+&~%H zFtUx{@~3P28Ss-jgg)mwVpIuvSn%PYGT7o5$fF0blA8XcY4hmkHENr6f?>Zklk0WpbApYy|gg z31#$6h^yOjO$ruEnzS!$c*~juo)!C#KnRkU|9wLNc|?HEevN#`%)!+pc>xo~pwzPu zen@Yo8Dw3=?*uwF!Rj>&!c41qoNW5;>%o9v#Qxdr0rPmn1-5ko#p$|m-G0nWM6--a z5y4}Zk+k~OW6AR3s-7w@sk3-yTf>jjgqF7Mg2XYiv?RqP{_~6Ph*;Y1hIPvnn}0=E zNtq9FV_p!0JTl(}TEKqcmJ&0m8vGW%L1$?qQOR}c%|FC8YM(FPup|}4-VK~N?+miz zY#sAa`oYowGe(2kv+T?OT}Js6pL@f4Wx58ZNrOG_`QD~U7wz@;{3kXg>jR;{nYngc zP$C%5mdE6X2YAwK;oES1`2g5(9}q#4{(K|IAHaUb=2`iL zQBbp%@eM~mqF%GTlUx#raiy%oo9IB9^S|NwrHk@0`#T^X5pW`LQJsqoh)e_<8xC z@AW?yvi_?N!kH78{qF)Dxf#Ql& z(m~2+S;_F@g_Oa8f6MCrR;7Y3SYt*y5FZj&z-7jLv+t3ywcbD->4Z#V-=j)H#Ln9- z>0DV^NJljEd!`L_%JJC+C98l#nz1(?w846LWd%*kb=sLfqD1+s@1W&9;%I0oESwbd zt)tDx!WBRB0=KC&7|@z%gi#&=U-+uu-m9!s3O8tMmo(v|nyw(#h7nIwq%4W$aQndd zSr5>kR5Sp@(ElZ=-DZ+dZ+c+>VoC#PcSfg1U}J zIb@RSev1#}NmF2YYjbgcY|3fUPt$}?a}>ojt2`9%nWz`Ctl=*7Nrd9dqNxC8V=wZK z990xR!i#WmPx5WQgx%7X+Vph|GA;_2Qi=y0?y(H4D|pFGB^62`!Eg@Cq>H`8w<7HJ zUCk+e@?9bjzB}y6(&2yc0;No3BED*qL?dK7(i>+B<)K-Pscr2qi(xV-b@O^ks&{>^ z34PdEoV<=dGPUE;!ntO7W(OWGn5|?_I#XpWrpy}E%P^MG6DA902k?)@Y*c7By6t_R zocq${!dB}qCvFG&gdZD%PP~0-9JSVdf#T0i>YwYMl(EggxBD=1ijanBCCXNHF!I># zd8|zc_gszj$cHtjt%F%VN+-@@fuew5FhJE_SZBK$r4rDXJ*Y zapykAXmSS21!>C2mS&5u-U#r%?@u8Q^_VSu#1*zT1Ue3RXxDqr*Dg&d7^!fKDiSZN zW`yPYyS9tuup=tmvnK zG!b4;IYXJrUi|K}wdY(kR+R42aPo9dylYv)-j1{#rS97pj#aA)K{|cQ!{1Ey8pw9)Cgb7K%doR_EP-SG@-<~&Z2L($S!G$y$WrS1g>px+;N`o z-U>Q{7{OhvLz47Efpx(*nR@j-6})qk9aMNsvF6R5sJyN@5j?|C3G(VF9u4q#`MGpxFnZ(2pN|k!KX~W{ zAAC0~XHFAd_8}?Lgae&I1%a*eJDTTAU2BhElMbkwo|SM8;Y`^i}}`3{x?f^^{;?8Xxn z5XWrO_}w&@-Wny-S9rQkNbTu%S3pkOK;5Wr6v)-8~T+m zP}Z|Mp@*(xFF-}it^z$E--*pdBczqJDKq@cK5>khY>n-sjujP&W1S;|^Li8dSdO54 zI##=wnNV+2laTL5V5K`aEQ*Gb8U0!K>=HCWyt$-69{`l?G{HBqAiivBOkdFi zR^)tUFJm|J=3%a)fcd5n@_;&Pn2S_l+UDs&_BPdnH<0S8ZTFd8+f4MS8(vgdu+|S8 zt_FR25Kyu;mVKFU>0TLBTG#MAJ{+K7A_KK_)UE| zL;7^~uXqg1vBnm(8GKcq8Wfe6tQOJEP}R0L-zvg?70ev^Wtt$^DtHa(^{n4~6Y9yq zHLN@1QZTLELbINM?EPq(++l>IKGjIkv!&gF?JV(y^5>*a>GqT%z-~x&qXsJl`j9XK zr^M_xY;$ zFvE7aso;EgrK6$18d|mKx-Qs9yy)6vD5(9QSA79!QYz>-uge@TP0gpEPl&pUIV$e%v!{ypjCX|+(1~Uq;>W?XyqYxIZm{N{_mG#ydu{dP_N?8! zWV;BIQ^0bde%{dRt26}70k{8!i)``mxD{m~P8cr|tq~kIqljDTH;x{exC~gMX{|=k zI}W3ioJU}?AAu3Mup9N30v>LlD1iN|m#&qU^d6(mWpWPLMRZop7ezHjk*Wz1%L>7z zRF*C`95WGsmZNhv^eL4h#IU|XHW)HZFz0GjS#gI1gVDBG2Y?|O2I;ub$Z*}(MHje_ z0(0@l8e}m8a3;Rz7`_jHtPSCVMmRn@Mv*x_1qyI$`g|1a3`P6FMAE6Y- zXZF2!9^Wgiq>-pKn_{@V)dcuEQ0a4okdW;0%ZKDY!gX#`tL%8IZziDx&}p&k$}M(l z#+-^6#yb9okD?}75iDIS(|OgCJWbFGXv8;ktSmTIrs25hMQGZG2t#`)ga3CRQ)~tH zZozK8h!28(OiCY3773)waRxzvITcB#nsWSpT;unVWs1@Rp+@ooKet%Tz)b0hYnfza zzK5~i4yoByhiBg?Demj&?D_Wo{eo|E#zA=nlGx7$VDT5}5qzN+Gg>VYj8RZ?z45pc zXA0}YQEO%4QlpB+7qXuzuGP%()7v!^;{86%%ou<3NdD{gwf2;D@82mbWyA^?ZwTHR zE1bZzotpsTSXc1lLaFIEsJ%C`7j|0o%|1SWtwyOqRDINw@4ex0plxKJ)PFK4xFmNn zjrG$SdnoTOkQ_ue@&`~yjLix(=)~1;$b%m2xHxD!CWs(!Sfk22QgCX{cG&P2bNE{V z4*wns``mF25t1~h3Kr5FuXVGPjR!0ps<-oM%F#A&SVEU@xVaQ z0uGt=ZO5c~poe~+XWxjibrtBdw2&CE%Oc<-MBFZ=tSmL|&D$(EJUH?${nh61sHkw; zGCkLn@VO;TO>Mt^3@XnMzYmdJGdKS3{j;XNYi9}bKM^0@^sV5X4NG3pFlf{fW)|6K z@IK2Y>z$kN*+aqiMs1!I-lb_=O}jPh7Q~&|J(zt@FyJ!tcJl2jIfI^qukJ6dIz59p z)jQq~EVy2Ak0-qDNl{!6BlC-#5pE3YfP6EP(wMOknW}F!VtaUrEuupq%@Kp`84J-H z+5N)-namKk9;h_%{g`FCsYma`v8Md&L*hj@1A6pM`bV=v0e_+hEU@`2ol>40+)KUCGUEtq7_bKmKmjjT#(1B~lpI%`W?r|0vh9F5Az*J9lQhJg2|^j9DOwm?3g{ zUuh4$YL_ihTo2{+sms1oq(MhD=JmTUkd4A#P~1vzos#*zam;6427BBDs(hnAi&Om? zQ$yBqxY^j2)#wh=#{@8w6Ff~@l9+Z5o>N6aTIx!--9vd|=OoNRXrVTIpN4I<-7L{K z-q7$x`O%M#-uB%-&{=EO-8Xep-X~p)|o9X6()f$lK)nfUEYD z*Jgltn1W3>!wLo8)7f44I>jCKdV18=Z991nY78;E2Q}qYNJ~8}{)Aq*p0NQ6i(+l= z1^dVn6UW2%%a>jeoW*wD9W`E5Q&a>TsL+kThIgnuyO=eSvADaqjr3i|L`h!A!LDPU zlJ%TkjGJ2-?RM>+G(591RAI5TGq>mM#I7!j{P&BB9~2dzoiKm@ey_=^Hx2;YcC@3` ziMYmVqvo;nGX8 z4#`-Cg*5D!mbj|{?Ifnn=UhXN-l}~>EEm|^D@yT?DA3YohDyU5D%IbVKj(Q*#0q9p zXS=VQJAdWt-pNFDUJ~`2@|~UX?*^yDalC}ExumDKp}THf(ZXnW)mQ74Yk)hILPb9T zgA$&nKcuFGXPYIQty#;AbXR0yd)G93UEY%-p_YY!g1z8Z&7UM9Yiz~1S+&#k@=pXx zW$~X9*xz4-f2LqY154Jxbw{w{g-N#Hh?cVR@kqeZ5strh*yuFDxU-)#WKa3k=xbai zO$abU2l)~TzY|ycC|}=n-_L6|;ygC7u0Xz%`MQ{TR@adqY5Y z*|DIafL|9ae@yy<15f`8ya3R*3r50I;`2ib(3Rn5CX|DQn&G1M+Po9pT2SQrvO+28 zw3CYU_HFHPZ+8SgscoM}#2$T(9948A(0L5#2pT|jn0Zlg=@{dajED5r<%KtLi|J>% zt%H^lB(kAwa1mlB*A%k$L*((*7h}W9k0jqx`m`1_r2$oF3FZcDGM93=D2V_r`46t=`$3v@E*n@lD!MoVx_#xNs4nLk3Z3cRoCQw|ybNTfv{Kd(5 zZU`5>=07xUK&yvJzme-+Pt9&LRIacQoGNJr&6obsP@^-h_M(89&?rFmfU4rvJaN|X zicU5t#1+5g+FBo2d*Hx3fj)L9>)n?fwwD8M`Ps0!xA$MZEL-7z>~dE$`{Z+ymk?43 z;A&wZ!h(bh{W5FMrP7H+TWcYY>~`eskcZSM&->=^oG9}WLsU=t7^JEb+DK)s0i_;} zhc$6srsryMJZ`|G)PoK`JKg8;X6ki(U9L#6Z1EcV<=K9+kv-avsl`OWj{k0#UAV}J z6InVFm@ZcmR`tDZM~x1>~t5I#Fekv?>px0qHppvf9Yb5lAo>mSfE>XW`&|^ za{m5GceGt}^z>?wm|$jd2s1^R%1y6N6S9WYf$^i>_LMRhaiF71^KXFnY*##Dtrn7X z{P-!@+Jri^c;XzrIKR?&O7{XUcFdwGSxMBEpT*!5jdFFYUXl{VsjcCD_x7XvI-3WV zbgIG{It0oXlPWCKyaQ@HcXV*Cy8ik!Tlvj%41)RbO z(*^HAWcSBLNO(6lXO8-XwSz4+9= z0NwKb6F96M)IRfUR3>E|mBM7)851c7AX95@aOVcIF8nnZE6YABB*C|+XY)576|&;s57 z9#<_+|8|fui?zOTHqWnw9$%voxY2>v5pKAVw%{OBvh#uDK$MVD^%d?UpqrS3tM=(( z^FZrKr!*Gyt{QcrXf zoxs&~MD6o_VM!2oSqS<@YDQz`jFgs}Zm$V9z!{R}FBL7M!bED{LqC zr&WhH)hv2Zwfd%8j59xsGl$+kY``L1Od9Y!So>(SjR1ppH?*q^b2Y1oLO7wq+QvB= zTO=w_@jKv=Ah)jnvh;lNp^KM@x}|eIEs2!c8{PhIs}t~5S4@|JtBy)6P{0wfxfZgm z2qW2cb+D*F2W`dGoaI62MAONr?)yzklMEh|p}TUPR>iQ=~LzghiDjoITu|Ng1<4Tt&3gxH8KK z-HkZbo_KfwH4#oUz{=Q4B zzZJ$nS2$zyPL5sL3YF4Z##csJGR{=s7a)rVp_ycFuw}~Zrcyx}9n*1?9YObsjTYN^ ziLy_zHd(iD5=&E7si#vVBCx9yCLk$#nOgmbGHgIks`}B^1o?s;fp0R4uSWy77_0uB zZh;^G*VS>0f@0+>;W|1LNKT;6d13P0R!3}M4Ej^`Sxd9~vo_sv^pTd;)fl6W%^h84S7L^JXb<)l7H!HT^C!?X?hQYA;IHTeRpmpbZ?H5P?SDcCB+Sqo?NgyI1G^RoB^%yGQ%PR~@j*y>E*=J*Ej(SS2^!mg8<+U4 zEM6-09*FVuOrJVW`~r?XinRbI@T{@nU^*ta+GgghpUV|hyZ|0hsZQ#k3h4|(7)t@T zIg&{p9mQuSDALDU@mIW(tM~5Ycsfz%x%w^RF&{RB-C23A%ts&?UVJ zlV%%s$YW5%?)dBE_%}$dyDSkTB{VY53k=u`&)U9`3YxfQ@YNrICI@((s&kT6vhUWY zG`VAVjjz%a7oW8gC>jGDPyj8Oj_iZ^UM+Z_3x{OW}Hh1No6Y$;E5Mb+2qT&%SkDqxDhouRi61j)n zd$RfZl8zcx*7E1lQ!+z{OL5UGAc_s~$HndpJr#>?K+BM?Tj;>9p zU(F1`RZYbjLF=z9k%F_WA=1lDdUK}WZ%#sEt=#aa%UIcdc~;J`Id z>3Ok}$&t?{qT_3=%lyGgsC)~)?E~>oUaXRSw^Pg^Cf*wr*|)cSk4MeoW_)HiKg^Y4v^Fc!R*+xEej*FU{P5OR;kQ5tJ}h zwuh)#4RXf6%8dUued2c?WS!fJDtXoA@g2?0gP_?u=i2xCWLHa(w8X@-Z`bc%Sa!;z zbv5Cel|K=`&U6;ccNZI9lZJ!J2@s{doLtTn!@v~`@NZe0I!r2*54PMv1uZrpw&u=k z*?KT4WoK;oy>?3?ep`U}1*L!{`Ni#Z?$woW8P$vhcSC|ztZhA^$EFEyC{OavuOsyBYu(eDdt(`AityxyuyA*F*KJ*6%2 zsIg=y&Yyo6PKbgKqF4 zAyAI{0&Q9PioFj?s=jklvafdDL7vi4#Dl^(tj6w#ADWI$$oHWZV}-*MLtKx^Xihm9 z2y&(XP6PX!)R%tXvn3)A&J7=6Zfs74&)Ax>1s~-%U%h?}dfr|CdJEY8z<>|q3(4_x zWD6Jl22R{ZD>`U879}0E@630MINm>MY>kGA$x81!>$(yX16><+HLKUy8)rttXOju= zBQp5#m%0;lGANF59~7LjM0M5aYWQQQC_doUK4niC2aBkji$;*^ z0gS+7Bm`f7^BhcOxA2M{`rm{m83q238fKQuObx!Y-F$M~aaYCsCvDBNlLYyP{{i^` zr6jk~zkcF^8tg4Vms*VXN@E5;4Of(8*IMW>$@No@`#*pmxGh`LW4I52enW z;RUnK_Wv^HC^=LTs(_e-gN;{}&kfq8@5rRQ8I&AIz)7O0jh?av3qID8KC1u$uZ;be`@v7}a7vNJtD`CQ;d+^5sfprYD1tqpka&Jjv z$TlG-3509g)%TndaJzyaR!Fd9Y0NR{vox{7W_0*?Hb`&ooco;sZs6)M-5sUvake{! z1=EVHucQF-ZxcZPSE3n`*r`Ek4%*601m z_QoRZ`?uWS0v>iX7n`u)gS`H`%=*?XV; z?DyHv-so#g)qRdB%cHTWeoy$wImg}ZBeAY0l{%Xyyg*fJ`cYKMn;9AmBqt~dJ^kNE zZBc8)Va^IPLl>Z}sFnT*tSwUsHZmu>tfX;+=hUWt=N~e@W;mw}d?wXAjE=C_OWR^bT-IP4Q`jr;B2}H{_%;CdKAo zVeLnTRW}W&Pds_x{YmMjH0XgO4Vt#mZqZma=W)XO$NPH{3&`w;yl+=e_>_eT6*7;`uy z#O5eW{`SpTz4WxEFiTYJ8P`8yJyjm}y`S;IA7~+e5B!9m*7UCz1CM=$tbv0?f*mcS zCFph_kF_&gP4+bR{M>RCVA%)(+boX=aa&uhs={`hAaEg$>=-SRgfbX;Gh#o z12{AfNv~!?lz}>l1A32O$?~nDOYstXsE?W5Ubwei>)>NOrz-4BndDWw*ifzQx#8i_ zbqmg2b=w!&S<CSr;nT?+ZJ1L~1#uNq2O=XCbyAjgIi#mRa13J_K* zawtwnoZxeJqk-|I0%mkT%dlEZzUJ7mMZPyDVcd~*_{wzs3yk;RK?jE^Jx~X(K0sgP z19J`{^T3>!kxIHU{}h9+;xE(Y`dIq|?VAc3xH(%gPhbr+wh`+b;-mW%3_Tk63DPOVLCKqXd!$6E1L`9e zkeM$17+y_l{yKwMm7zSo_UWjQcQ000lHxwfI4f$$h7PX?9YcEi{2SVDe}p^rpHnby z&I!dG1IBYD>x|hC?e@dk%$hUkW|*KR9_;6f%!xjCMjXY?>!1oEY+szQ8WsJjPBM4; zG({*Kl${hGu15z(gb$3{!#vR5Aq1;`w172y18!~0Bjcvt4>_GA$KjgbVlVMY@cgRtSN(Bugf~m-oHe7 zSEMo1=AQj!aE`jM;zD7-#x$t${(mJ1eoKrmFi<~2^3ZZboV)KQNE_}sRx}&xogQ6* z%ia#ID%M72F8fWLW{p)F%k9PsF0_a~*bV4`v|AYT5|~5-`g8kQ5Im!XYiH3{Suq-) z9}PL27|yoE*68jXt1dsz53@Jx-A6V`zAI9FvELcFdzav+59ZSxf2(fL zI0$Bf?Q?}aGbV;>VoMOo+mSD48PpR9Ra4Pvy*G_bLb49mm#%k#8}km=jd5`xNR1JG z&l-O0aZGCxpZz+23litq0?^qbm9AW!l}f@qLx*bw7${qgd8wD2eF({iwbjEMniDAU z8+!T5G7G(j;qQ8o;!{+Gyj7xww{m6TiRv)7U2KiZtzD&QqUkh^zYaoQt0i|Hx-69) zY-8;1k;e$@VL`53|A*I5><=2^H+pE!cgxe6t(w@(N@hi8sr}GFb8NTdb(wI|jq0J< z#I2V5hS$U-Kl<**-QB-a80c~>rh53q!nei!cY+cRE&g)h_vaf3S;dK-!gJ3^` zd`%lUaH@!4qThoo&gENm41~3J7y7YWr6ZG6B*k*Hr#JSUXMr89kJ6fx!i&}Y#WS8} znB+43ws>@Ns-^%n@zDh^Q?zMFOr*vsLlpA`T-FDjPa>*;BZCmemnlYT>qqZQe$A*v znm;bt(lATKr(PrvthSvCU04yj)i5(X#y*0oq{s(YEzWQWnDwfn@v~;&3?uekg~d6( zRe#L_WX4y668rSfKZJ!{rFIPfS2}rykU4Ub$nH%zt4F|%0PQNA{wVycHuc6rbE(Ck zO}8^0$n@s?G=9dKN}U5bO0eaeQhZCv)!4wM)_Q3}nVIC2{Kyr;5`Z2fDyZvtoG1X>_+k(?nPnDU zIXX8~-Y~vhjQDbHUJ_es&M-*$GCvvHn3;K^g=k%nS|n$L01a8GD3Qe|4>C0C+{AHL zVp{pR-(%O^|8pw)oly5J?9Y63vkv@t4dM-=`HbOWoMa>OX#i;xlB0Ax?rXy`Lz)@a z&5!Gp3O>8bEgc;gMu(&0YY?aIi+S=uo}Sd4con@Dd?qBNPi2ntlr(~9aEzNI{0lcB zL$<4b?Sgc>G}0-j^ga07;`jpr1S|AaW;6wn!Ufdfi&3Ge0l4jmw@-U3Jl($Xz3fo3 zFF%1}!Yuc<2XHBj$i*Oq9n^+Gb-Kns$Oo$|T)wFzvpmOU{5LO>Opx6gVvloH@ste% z*JUNc+8zo2F0brwxK&VJWW2#l(zd`<{kG>?KX|XYU1>|dOCY#nXVugWncOPZY5y~E zSlmm$qx;4ZIOjck_R->52Mk}cttOe2N9#lDtD+vAV;jOfD|U~+5fq%p7XFqCJhb=w zj|p3bSJYof8^yVyDQyn-B3gxQ%|9G?-#y+)S^JOl3k`eP%q(|Ax-|VAcR}p0!?HWs z9~OBFY}h@V#^}7*N`mB$hPpftCs9j1oPQl&Dm52J>4VS{8_*I~8*{YlXcS;u} zp^K2owkXvyJ_8#@_rbp5#g({ReGLzBgd6n{zHweL=)a6*XDE?KW!OX=sEZaA}Gra_0&k61~8r zf!8Z+8R2$WkMDWAeTIu{{DORqN7_L}s(R5^##nztce4EF>Ko6c8N98yL=lyw}{wR;(fHV%NKAw%E>Uid|lB^Bznom9t3Hc z&({@k%#oFogU_RHP|Wup8?iyRxlG4n2!r|mwHx#M%KJ}fsei|dagTH-XBKcR^#t`+ zmCK3^E76cEjTG=GAC5wg_s;Qc^(IKP%E@}k*v?Ov#9JkpNwe$j5fVOM+OlECvGfk( zQUpRx2q352fo%uMNN-y31fNYQG4pX2jge dz4gjlJC1J<{2A%??|=8F^xXgBlN(-4{1baNi3