From 06dcf8fb4f32d5b0f807440e381586dccb91faa8 Mon Sep 17 00:00:00 2001 From: gdamms Date: Wed, 11 Jan 2023 13:18:31 +0100 Subject: [PATCH] feat: draft fvt and ilv --- .gitignore | 2 + fvi.py | 141 +++++++++++++++++++++++++++++++++++++++++ intersec_line_voxel.py | 102 +++++++++++++++++++++++++++++ main.py | 55 ++++++++++++---- torus.blend | Bin 2141928 -> 2126784 bytes 5 files changed, 286 insertions(+), 14 deletions(-) create mode 100644 .gitignore create mode 100644 fvi.py create mode 100644 intersec_line_voxel.py diff --git a/.gitignore b/.gitignore new file mode 100644 index 0000000..6c99216 --- /dev/null +++ b/.gitignore @@ -0,0 +1,2 @@ +__pycache__ +data \ No newline at end of file diff --git a/fvi.py b/fvi.py new file mode 100644 index 0000000..6a1e26f --- /dev/null +++ b/fvi.py @@ -0,0 +1,141 @@ +import numpy as np + + +def fast_voxel_intersect(start, end, origin, step) -> tuple[list, list]: + """Compute the voxels intersected by a line segment. + + Args: + start (array-like): start point of line segment + end (array-like): end point of line segment + origin (array-like): origin of voxel grid + step (array-like): step size of voxel grid + + Returns: + list: list of intersection points + list: list of intersected voxels + """ + + # Convert to numpy arrays + start = np.asarray(start) + end = np.asarray(end) + origin = np.asarray(origin) + step = np.asarray(step) + + # Translate line segment to voxel grid + start = start - origin + end = end - origin + + # Initialize list of intersected voxels + intersections = [] + voxels = [] + + # Compute direction of line segment + direction = end - start + global_distance = np.linalg.norm(direction, axis=0) + if global_distance == 0: + return intersections + direction = direction / global_distance + + # Compute the sign of the direction + direction_signs = np.sign(direction) + is_positive = direction_signs > 0 + is_negative = direction_signs < 0 + + # Initialize current position to start + position = start.copy() + + # Main loop + while True: + # Compute the distance to the next boundaries + next_boundaries = np.divide(position + step * direction_signs, step) + distances = (is_positive * np.floor(next_boundaries) + + is_negative * np.ceil(next_boundaries)) * step - position + + # Determine the nearest boundary to be reached + boundary_distances = np.abs(distances / direction) + clothest_boundary = np.argmin(boundary_distances) + clothest_boundary_distance = boundary_distances[clothest_boundary] + + # Check if we are done + distance_to_end = abs((end[0] - position[0]) / direction[0]) + if clothest_boundary_distance > distance_to_end: + break + + # Update position + position = position + clothest_boundary_distance * direction + + # Correct position to be on boundary + position[clothest_boundary] = round( + position[clothest_boundary] / step[clothest_boundary]) * step[clothest_boundary] + + # Get corresponding voxel + boundary_reached_negativly = np.zeros(start.shape, dtype=bool) + boundary_reached_negativly[clothest_boundary] = is_negative[clothest_boundary] + voxel = np.floor(position) - boundary_reached_negativly * step + + # Add voxel to list + intersections.append(position + origin) + voxels.append(voxel + origin) + + return intersections, voxels + + +if __name__ == '__main__': + import matplotlib.pyplot as plt + + def update_figure(): + positions, voxels = fast_voxel_intersect(start, end, origin, step) + + plt.clf() + + # Plot hitted voxels + for voxel in voxels: + plt.fill([voxel[0], voxel[0] + step[0], voxel[0] + step[0], voxel[0]], + [voxel[1], voxel[1], voxel[1] + step[1], voxel[1] + step[1]], + color='#e25', alpha=0.5) + + # Plot line segment + plt.plot([start[0], end[0]], [start[1], end[1]], 'k-') + plt.plot(start[0], start[1], 'go') + plt.plot(end[0], end[1], 'ro') + + # Plot intersection points + for pos in positions: + plt.plot(pos[0], pos[1], 'bo') + + # Plot voxel grid + plt.axis('equal') + plt.xlim((-10, 10)) + plt.ylim((-10, 10)) + xmin, xmax = plt.xlim() + ymin, ymax = plt.ylim() + plt.xticks(np.arange(xmin + origin[0], + xmax + origin[0] + step[0], step[0])) + plt.yticks(np.arange(ymin + origin[1], + ymax + origin[1] + step[1], step[1])) + plt.grid() + plt.draw() + + def onclick(event): + global start, end + # if event.button == 1: + # start = np.array([event.xdata, event.ydata]) + # elif event.button == 3: + # end = np.array([event.xdata, event.ydata]) + start = np.random.rand(2) * 10 - 5 + end = np.random.rand(2) * 10 - 5 + update_figure() + + # Define voxel grid + origin = np.array([.1, -.3]) + step = np.array([1.0, 1.0]) + + # Define segment + start = np.random.rand(2) * 10 - 5 + end = np.random.rand(2) * 10 - 5 + + # Plot + fig = plt.figure() + fig.canvas.mpl_connect('button_press_event', onclick) + update_figure() + plt.show() diff --git a/intersec_line_voxel.py b/intersec_line_voxel.py new file mode 100644 index 0000000..ff3131e --- /dev/null +++ b/intersec_line_voxel.py @@ -0,0 +1,102 @@ +import numpy as np +import matplotlib.pyplot as plt +from itertools import product + + +def check_line_voxel( + px, py, pz, + dx, dy, dz, + vx, vy, vz, + c +): + """Check if a line intersects a voxel.""" + + # Compute the intersection bounds + kx1 = (px - dx) / vx + ky1 = (py - dy) / vy + kz1 = (pz - dz) / vz + kx2 = (px - dx + c) / vx + ky2 = (py - dy + c) / vy + kz2 = (pz - dz + c) / vz + + # Order the bounds + kxmin = np.min(np.concatenate([ + kx1[:, np.newaxis], + kx2[:, np.newaxis] + ], axis=1), axis=1) + kymin = np.min(np.concatenate([ + ky1[:, np.newaxis], + ky2[:, np.newaxis] + ], axis=1), axis=1) + kzmin = np.min(np.concatenate([ + kz1[:, np.newaxis], + kz2[:, np.newaxis] + ], axis=1), axis=1) + kxmax = np.max(np.concatenate([ + kx1[:, np.newaxis], + kx2[:, np.newaxis] + ], axis=1), axis=1) + kymax = np.max(np.concatenate([ + ky1[:, np.newaxis], + ky2[:, np.newaxis] + ], axis=1), axis=1) + kzmax = np.max(np.concatenate([ + kz1[:, np.newaxis], + kz2[:, np.newaxis] + ], axis=1), axis=1) + + # Check if the bounds overlap + kmax = np.min(np.concatenate([ + kxmax[:, np.newaxis], + kymax[:, np.newaxis], + kzmax[:, np.newaxis] + ], axis=1), axis=1) + kmin = np.max(np.concatenate([ + kxmin[:, np.newaxis], + kymin[:, np.newaxis], + kzmin[:, np.newaxis] + ], axis=1), axis=1) + return kmin <= kmax + +c = 1.0 +points = np.array([[x, y, z] for x, y, z in product( + np.arange(-5.0, 4.0, c), + np.arange(-5.0, 4.0, c), + np.arange(-5.0, 4.0, c)) +]) +while True: + + fig = plt.figure() + ax = plt.axes(projection='3d') + + d = np.random.rand(3) * 1 - 0.5 + v = np.random.rand(3) * 1 - 0.5 + + px, py, pz = points[:, 0], points[:, 1], points[:, 2] + dx, dy, dz = d + vx, vy, vz = v + + bool_vect = check_line_voxel(px, py, pz, dx, dy, dz, vx, vy, vz, c) + + # plot cube + for i, (px, py, pz) in enumerate(points): + if not bool_vect[i]: + continue + + ax.plot([px, px+c], [py, py], [pz, pz], 'b') + ax.plot([px, px+c], [py, py], [pz+c, pz+c], 'b') + ax.plot([px, px+c], [py+c, py+c], [pz, pz], 'b') + ax.plot([px, px+c], [py+c, py+c], [pz+c, pz+c], 'b') + ax.plot([px, px], [py, py+c], [pz, pz], 'b') + ax.plot([px, px], [py, py+c], [pz+c, pz+c], 'b') + ax.plot([px+c, px+c], [py, py+c], [pz, pz], 'b') + ax.plot([px+c, px+c], [py, py+c], [pz+c, pz+c], 'b') + ax.plot([px, px], [py, py], [pz, pz+c], 'b') + ax.plot([px, px], [py+c, py+c], [pz, pz+c], 'b') + ax.plot([px+c, px+c], [py, py], [pz, pz+c], 'b') + ax.plot([px+c, px+c], [py+c, py+c], [pz, pz+c], 'b') + + # plot line + ax.plot([dx, dx+vx], [dy, dy+vy], [dz, dz+vz], 'g') + + plt.show() diff --git a/main.py b/main.py index d81b802..5b6761d 100644 --- a/main.py +++ b/main.py @@ -1,30 +1,57 @@ import cv2 import numpy as np +from itertools import product +from rich.progress import track from matrices_reader import * -VOXEL_SIZE = 1e-3 +VOXEL_SIZE = 2e-2 X_MIN, X_MAX = -2.0, 2.0 Y_MIN, Y_MAX = -2.0, 2.0 Z_MIN, Z_MAX = -2.0, 2.0 -# grid = [[[ -# 1 for z in np.arange(Z_MIN, Z_MAX, VOXEL_SIZE) -# ] for y in np.arange(Y_MIN, Y_MAX, VOXEL_SIZE) -# ] for x in np.arange(X_MIN, X_MAX, VOXEL_SIZE) -# ] - projection_matrices = matrices_reader('data/torus/matrices.txt') nb_frame = len(projection_matrices) -point = np.array([1.0, 0.0, 0.0, 1.0]) + +points = np.array([[x, y, z, 1.0] for x, y, z in product( + np.arange(X_MIN, X_MAX, VOXEL_SIZE), + np.arange(Y_MIN, Y_MAX, VOXEL_SIZE), + np.arange(Z_MIN, Z_MAX, VOXEL_SIZE))]) + +background = np.array([18, 18, 18]) for k in range(nb_frame): - - proj_mat = projection_matrices[k] - cam_point = proj_mat @ point - cam_point /= cam_point[2] - frame = cv2.imread(f'data/torus/torus{k+1:04}.png') - cv2.circle(frame, (int(cam_point[0]), int(cam_point[1])), 2, (0, 0, 255)) + proj_mat = projection_matrices[k] + + cam_points = proj_mat @ points.T + cam_points /= cam_points[2,:] + cam_points = np.round(cam_points, 0).astype(np.int32) + + visible = np.logical_and.reduce((0 <= cam_points[0,:], cam_points[0,:] < frame.shape[1], 0 <= cam_points[1,:], cam_points[1,:] < frame.shape[0])) + cam_points = cam_points[:,visible] + points = points[visible,:] + + solid = np.invert(((frame[cam_points[1,:],cam_points[0,:]] == background)).all(1)) + cam_points = cam_points[:,solid] + points = points[solid,:] + +for k in range(nb_frame): + frame = cv2.imread(f'data/torus/torus{k+1:04}.png') + # frame = cv2.cvtColor(frame, cv2.COLOR_BGR2GRAY) + # frame = 255 * (frame == 18).astype(np.uint8) + # frame = cv2.filter2D(frame, -1, np.ones((5, 5)) / 25) + # frame = 255 * (frame > 255/2).astype(np.uint8) + + proj_mat = projection_matrices[k] + + cam_points = proj_mat @ points.T + cam_points /= cam_points[2,:] + cam_points = np.round(cam_points, 0).astype(np.int32) + + for cam_point in cam_points.T: + cv2.circle(frame, (cam_point[0], cam_point[1]), 2, (255, 0, 0)) + + cv2.imshow('Frame', frame) cv2.waitKey(0) \ No newline at end of file diff --git a/torus.blend b/torus.blend index 0d1bfdaf6fd27e5884264230f2acf446b446ea29..173b5682978ab3f0b44fe89fe2655290ef1d4553 100644 GIT binary patch delta 7085 zcmcgQi(3>&^4r}rFf6NpEb?$!KtVwk6h%ekv4F-k0mbtXO-Mo_##JyE6_qF|Y98xF zWAN0%JB^x{$X!TIOmGkvd5FoHV50dvLwtNca;N7SiJCWeA^NMHnRWI4fSd2@>FRoQ zRdscD&A2abAg$Ndlh=(zR~=GJFE_4Vxnm`XB)T0tNObT`w`gP9(U!!B>I~J=veAX(XCx zm2WR2gz1QmDfF3tms5A^vCo(ZO{)u!lqbq|9JjQLkbq`ZYsU5e<{A269&P z7ZOgYwkDi}*1wU726yh*c8-b19-aKJqW8@^Hx;5#Yx}+p zy7Nee>@X4=F-SS8xO4zH3L{w@d-LU!oFNTbO+=J)29iW-ko4|j%Cy)k)54h$A_=xL zYddy6k1~({tSD7v1h9mobnGC~B{&}Wmvg6ybbt^^!qoU!{D~x6BFS~?cd# zCz4T67)`c7WDE)QGQz|d;udNzH92AAU{sSjJGM-fYED6~O);!nl&*1P+bCpA9n3x?In^nN@S&N7!?&Z6Cw~aya$$V}wZQ2p z76#%N8bzggxM{{g))=}Cx(73f)`Ze1A8E6~^JXaWd^Mhp79IJVVSN~v->T4>!srFj zvFwsl3&nfkT<#5p_C+{7C{Z>4B9x%7;(YZHT>PNI@=XM)j>$yT(pQ`oGm{wQ%Y$ei z^#~4A!_hB_Uh$yfx$~xFnB|A!s0f+Sw9`XMP`NW>yxhZ8Q~_Qruz0YSHBrSc9^z%S zaK$gi(yJbINgS~|&WlLalntbz&^eTz_SM{j`dGmNox{AGBf036czX37%YF;7KwpBF zWfW&ICVCO0RqaU@uSuz#XJV3+-N+{&45u>WD z!C;)~wZP2vTc#=WtRYfdwY47@qHwyRLV6lyo=jG8Mbo`9F^WyU^CAYTZfu$1b#(~m zIWW`96RYwx&+_ubA&-1+HXabLBmtdCeH@($ibE*B$Sdf}H91G|2s}YgL(VXk3TvM9 zYVlJQ?|6z{5%>ba{!{u9a$PYGz}Y00Dqntv_IYDrn2(c?k*MT}>8LaFS=xziMz0o4 z)z*QKv4F-yr5zL2W{0&!)C`*!&~A8fA+o+`$5-bf>L)iZq;@=0 zWnR<-?nUQfuUi=$apgHLB2yLm?e|`xEFK`^657NdD~t6X>h`;>6jyf_Gu=LQD=P0W zN2*#~Eib?x=yLj;o)>FUi^vc9Rm<@Lbmi_+cxrNMPr}G-W`ax0=_P&vhQOXjnB;#| zyBx4rD>T%Z%&$ma{{?1 zNa37w&@JxS0IdaTc+MA%1kd5@X#?F0 z=QiNd6t7mwavB!qQ8CKrF%PJd9ziiIqGuXh^RFpORQtj;Ex;T9Pc1U!?v1pFK*^(6 z)D^r$=jlam$M;G;G`&RAja}j8%Gy+XvPZeOlYWp+SAS#i>3GpwgQ)j@rSLjgqPcc>`ehzvNRn`-hz%`GSlCo ziP{!xqV9=U+*fVHDr?Iox>qY%`{c$p*5b>w1rG&NIXeT!*N}=Ugy+C(bfmJriZ#7x&Aur#M2_4_-=gSoTQ&Fi z2Sh9D)b*{p5$t;|xRhEyGk&RRBu@|i)K0pfB)mF82UCAgwuWas9)T3BarqMi)hlqF z_H3im*6r?F*D|l@w`0?SyN+?1BdH+0Vqa=Y^1?CMeB5|9&bj2-G1--Gq$*>MsY0(&NM9cXt&pLg7NVYEAo7#d#v~Z>TIlEGygGS`s zotDTiF_<5bb_^b`vlvPj^j%gl(lf}>RcEyne3!j zS3v&`YL&QWZ3|4p&^--@}AVZ}-|cCOf{9xp-;`q7hv zICl=FiH9J-p~s0JJ9AkE5d@ts;9y@Kqk4y?()J4^XEU=L`w+ugp=vDn+t_>P^-&uu z{bToAcGxFKez0Q-yXEgVtXk3gz7DoN%oHmO-ox%Y%zE!(SCWs`Sq3&K!nkLEFK{hUTZ>*H(yOf~rWK-K$H0O!-}kiKSHrjkHu zhn3Za?awGVx?~>as3%eUprC}sGTn|ma6ZGPKywKz#Y9M$e;>6^2#eP<-+RtUkh%dA z(zFq`N7zPIqNFxAUHR)F_iSWcP4cC3p+V+9H+-c(;;fEg3y$z_c>E@?h}HU0^aw9S zq<^=|3zEdY_YabWFHkVGTCaop%@}|4R^Q?9{8M6};Q07_gN$2Z41vp;;z5@7;|zXh zM#=AFiDn{5$9tUM$`*GU`7b=K+tWMj(~rQR2x5}ATSZA^XbJ3{BrfhBAFod7QhJEF zNzZ zJ>@?gCWZm3mW%!Tlhjp+$DUvG8L?8Ps)ip|Hg177-H%wg@Mw177QpOYR>DvgpXH-Ro9xv{#@R!do> z?p2~uP#`V?&Pq`N(`s>~pzNZ;)ncbX!hdjF5AT(I1rsa9IRYQPvaeFa=7Y@_n}p3z z-dCwFoCXzL`u_5rPxS+$$(OGL_2{wrV>4juhpj)h0BnKSg0Kz17L3h^Ed<*@Y@yh~ zu!Uobz%~e*Ne=3XT<}q?| zwMXHb26L7meSZV>9zp)zNG{2tn$NEM8R|G>`wpo+N6eY*HfNJHpiD-Z<2a(~II5X5 z$=paMLS>CebSsl1M5Vd#O(Ra^o{rxwmZQiXh}M)H<=W5ic{tRF$2s&MM5P%pyvdwJ zT(GdooGtXep*sw9O~_*3@KXcyG?|@(YXP757~LIY>)T@E;gUx6DP z<~2gi;PR7Db`q_04~D;+lKca`g{Hzc^%H$Nl;QRoG_vLE@*#AC|WZN)ZqoGm*dU zQ{<|hf?PeH@-#H^73}0X;s{a^@^j73Rla?vz;~VIan$`N3_pzv6>0B|=#Gbb{`cmO8mcOu_ElvFQlS=J(&b#jc`l-j~a?iWx-}jt* z?tSmw+5hh)EO=%qTWe(*H8GQTQpS?G8|JbEma$<2ONw5dd+?dXCX*`vXOGqw6d$Y) zum~mo#(_o}kWAw8zuUr&DRAalKHB2j ze(Og964-X@M{tMn9{!px>@Eg7mhloxAh1Z&RI!H^8%veF6GvWzmVZcj{wgagX0XS` zRN4FdKLg-?n)i@2u8q0;?yfAITYQ)P{Sg_b9-3ZfMHb&*g*)dKbf;Miq@PFRHa0t^ zs9;P%3C1JBH<{*flZnADg+*K3Jq`*WrsMU{WM<>#vX|D?liI*^HC6CX+Vm3}z#P(f zG$NT;Te3y&yiZ*eLD*X-->TP=Uw=WHrQnQ{+PxR=wIs^fKK=n&w?oBjvuY=2q(Bxa z;YLQC)Be3j%P3$*Q&WSV9qh7fvwwvdO+=JLbwpWr`ge?tG>mOeLSsS=jwi9pmO!LH>HU&WVqHRW zGONoph7Fj+{w90;< zJGoruz%blBF$|(Qh7mmZ+s=Y*&aU-DXF?6S_hMJ1=wIr==46L6kmzX2`_JDt<;5jN zavY92(UezHgP)!3dkfVHv1&(uifK+`Q$xfEJlva&Rm*Umv}&5^z+9b9W2)Jx>)DS2 zOXgWM?r^BEGI?V>lQYRkM?AsbnQp@=hwZ%)k5Ka;@69fAOR!edqW6p_Bqz-DFm+=; zR-I|E+NKFb$enQ+W$J9p+`@inUZ}zA{ji$qx?t+I_Qz=;sv}nQcji3!fRVcgGtih& zgQXc*O*WS>;lmlEns>rVstiO^s;1^2Fp%9i0uu+Zq5kQN+d{x*W0CNDD9@;IC(+b}!G%zM3Nps9 z)37a!_l8adSh{Dz_V27Ez*SEan)n*&^$4ClH2;of*|lUi^*GKD)}MLho-Z0>__ zAGPyDu%_@hc-f9)NBi<_P@ce(2(K#@y}zHp*TY96Sz3+TtoUkv{;nC8ICvIJ9l+z^ zGY4M>7q+OqYTOpZS4+%UcjF3pGtsd6T_XQjQ8y81SCXOjYZ9vM`cCOtIDTcJRN3?LY%Ii{HRRhaO{d4UE{Wcs~10K zNTvzN2c5jx80xJnr?n6D5knUyAc)Pq`4M3Ac|wi5D{Si1ZmMrhW*>Mjy?vmc(86pr zhw8xKl|(rmn)|g2Wsug5{`_2<$mBGZPVu&396%x%WZOWNFg!e~9pONd89%UHW{{Bi z^B~>}T%E8YwhiL-a17sK2zT{hjHIYDkCk59)d*YfX1ig_($GlC&X|$N8Xt+7-gY;O zgIjUZIxO9^C_WJ0EW~`yjuP{Jq?noEM-TfOtQm*H)m;pYehJb#sJMrzxN*}o$Q3lN zHPJY)XweQyBUvpQhIkdy4=>6b7Cnu>&!}Kpfn6FME*{X$G`P zJ2nD~=|We-z@P3z_Hd!Ofiy35!&T8YmWRrLs+@xqtpgT~Mxnbk8(Kjo$Dz#TK`7(v zoU6$s!c%d49z>S1B<9hqu}$#YVBQ}WykiJY^8afvKP5q_LloipTz(GD4%cq{vWxFE z!}9TtEdSNJ`Aw#gFg*=3>$-;@g-d@6je*f!aE-Y}pxP_zM3H0a+vYe?b4exK41BJm zOU;mE81+Y7qnq#LM~wU=)2wg0uU)IV(5fHZUKuGQbEKqQCRNCsD{YtQfinIZ_w&9o z9P5v7G)+gs!fV_yF#aJ5>RvVwW)Su@Wmx)ArD}YQ+XmI+kQY^gyzLJo# z16$2mVDa!+|Aq0qnZf3F!)WFEyJR~kwgOP$U?l$1yl{=vxmSq|7|@~L38VL3iO zv+a?f!ugr}wEjqVns0@oFL0|K{|pvf^)t9Ieb3-l{rwrfUw;^kW%w|d%|FAy>_epxCjL}<@K5o2K8}|>PPz)r$#?v-%+_b;3s<{5|mrCGB8M-c? z!ykpTaV&YmmyS*w=kjv#6sw-w0dn(PguHeRKE*c8q8Un(#lp67z|VK07#<^CXMO zq-Y-mxu#^k;g};a!lcg8Ki05^weiX(k>sQmJXL#m&5(R>dIh(Fw~!6QM-V-BUn#^T z;97~6chpO^jWd`G-^{~*H&${dod1_}46gYwSe$2kkXP{x^5p-eif`wzqZ%i)u39{M zdn3*#C3WBxc>6ux4;HTBqZjAb{HI}rG`t~k`ow(7kt;AI|MB^+@9li>+)%pDf8o(R zpO2oOzoOq*?ViTydg;X2gCyBIKNv~_k9g)=)S@k6ZavAeVh8dod4wG+Zyr?oohf44OO$Y?1`E- z%e|%l*aigAK+fsPJur3{f*2Oov#w!hYf64;d5=96TXWyd%^$V*=7&4?yl`mA(5!nm z3|!e1b@;U{8yiR_ zjCo(>u(W+-S~xP99R9Ow`MEh!%oJh_VXr33kk5EuzpP45v*p=cri-*ByziJULQ^Qm z2Ui#AGz3&`x2Qj#(B4L=7oVg-{1hc=^%dij)Dt$L5i%-y9zNQlF?yq$$uFsHjpr_u5U-N^De zVcGY+_*@RG5tcpg8kU{&{}~|5o(1Addy=1C))=gb=C)VgjFlK{x8jS=;W+6rM~LUY zOE7#GkRu%r556r@k=4v;_Dt1U-a1PU-S!gXp0O@_3!gV(9kkcR{Wduq#+Dy|I97^=Alw6vv{tTY189t zEyL-V$rJoafz-_SehGdUCcTMSd7)5x!lEt~hMg6LVb3#bbi*@?q^serZ;9A^@07t) z!$oWf;nzeCRqXZoMSX2nmKb)tnjH#lc(HT}(mj%zbT6xKD=4)^(EzGH92B*Jy|Cq7 zs$*J)wzLzE*#@NjxHtQ6dN5YpjdDgveR?G0PD`z1=OJx}rkyylOZz?T#L51u5z?v1 z#p-80-Y1Z$2!i3?{u8&~PrBp2=^wZAciwS7>GRw98F%Pc+%CV^)N=cP$#jRp%scp7 z@8Fy6aQMv~{E|EL-M8^$nW^Fq14-h3Lz5%MOc@cM!|?ie=^abKD&CLww8j1pj7wmb+~9_IsQpg zrAi)hi|UHO_kWU(n}Z0LP@SxWj~|nMMjHpbxT@5Wj3~yBP$(G%#LUL?gTP*Z@VzWRc($~GsM0|| zkA)h1wushviS-i#8^B)McBH_nR!Ux;|MVj1`~)QQFP(FF{JTK$B}Fk~o@V`$lk3c7%yVC32GWWd%GzB$Zd zBsu|Vc}b3m)-??Y4MSukv&u9_8y8`{Rje|ri;aYKbWMpj^o&@jPtXt2DRm+EDj*Jh zg@l^XaE0`YOWhx)?X5=dv5(R!gDOC?|3~dOIL71?gsE)gFzpDv6?=%xl@bQ)8d-XF zRfo_tJcvH4Nc2qwT3=BNqoY{gZU)p?8g&<=P@n*Ncw{__4|Pa!h&>LXY%DR1B#)zu zA7C=RW~LL5c0L>-8O7#4lW!Xo8Hdz+4+ zaE#zJt0`c4N(z$(bg&Z^ud5kmj7cH$WTpoIXNEYw?rXF)jIahjiTdZ3UB0 zvm!1v6FuWk8DlK1Q^l_cfk`K^&T-d_uEVr^Xe?BiJe{o2)ASszq@cy=PeQL4S#&q$;xJdTk0z zf=S;>R{Y4x4zEKrB<=ac1Q*Rpk7#2c5}{Je#OM{StdWs;ADU8G3=f##vv6evN3;)F zlv_;dbZL*jtdmj|7p({4Mom(V#Ct7Pe^Q!KW{Eb!+N~Q`F}&4-l+CF~XfO;0^0muQ zczW7nUMTJ3aDi{IQX--~Jy;npdmqTrJfy;;yObo^cl&K!+ieSR1N!X%_*X&b0HNCy zrKAfrjnXtqtV&;rWO2wDN^gCbei)LYnM)~>X#}>q6st_Kup7rD3>Gv zVJ|pHLH1abCBo_g`NVJ}e1!i?cRK=QP{d*i^gWf7K5!Twr<#L*kZB8 zVT;FR_g~+cFz#b~JuRK<=w{Cb@k^1>--BUM!^fElm2(~CviIV<^&mg%aLLXqmiM8E zsAsp}&nlJABG*^`*G7=%;aDw;KrT^b-DV$#%6Z83L?b27M@qegx|Q=CE*|v5!TFA% zvgg1!)Lek_C#R>PCQ)ijf0QCBknL%N#s!W%>HY^%tQ-l;*uzj%PEvnI%}N>~6OHBM z=&R?Exe%GQ!WYmvQatz&99)Qw)N%NzmlvV<+opG*Xb}o!Un4X`1tKOLf<~lxuoU=W z$51@oDp>5ump!h~51?`}i9CS61UD{r?2&!sud;&1t{~7xqUtUqB`2RB#MVXk|6}JxR~8Mxx68 z57EkNNOR{QC|rY@zEp-YxCYZ7EWQQOb2#KlC&bU8lzQ_lEW;tm<39w4i07MQ!$SKD zMYF8{rHE2rMXK^IXu@~S0@7NPuZ^A!g-FN&%2efAl=4g~#2HzO(R+80W({dxK6?Nv zkqYxo5kC%(fKO^1ck#ej;OlT=0$Mq*LnU=U^8sjFM@l(J$?GYI{>}cv^^SdU@{Fe^ Rr>zWm%i#}sEAhj_{|{(~Ql|g_