{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"pygments_lexer":"ipython3","nbconvert_exporter":"python","version":"3.6.4","file_extension":".py","codemirror_mode":{"name":"ipython","version":3},"name":"python","mimetype":"text/x-python"}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"code","source":"%matplotlib inline\nimport matplotlib.pyplot as plt\nfrom matplotlib.patches import Arc\n\nimport os\nimport math\nimport numpy as np\nimport pandas as pd\nimport cv2\nimport random\nimport seaborn as sns\n\nfrom skimage import measure\nfrom PIL import Image as Img\nfrom IPython.display import Image\nfrom math import ceil\nfrom joblib import Parallel, delayed\nfrom pathlib import Path\n\nfrom sklearn.metrics.pairwise import cosine_similarity\nfrom pprint import pprint\nfrom scipy.spatial import distance_matrix\nimport random\n\nfor dirname, _, filenames in os.walk('/kaggle/input'):\n    print(dirname)\n    continue\n    for filename in filenames:\n        print(os.path.join(dirname, filename))\n","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","execution":{"iopub.status.busy":"2022-08-16T12:38:09.911787Z","iopub.execute_input":"2022-08-16T12:38:09.912197Z","iopub.status.idle":"2022-08-16T12:38:36.686773Z","shell.execute_reply.started":"2022-08-16T12:38:09.912147Z","shell.execute_reply":"2022-08-16T12:38:36.686026Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"ROOT = Path(\"../input/ultra-mnist\")\n\nTRAIN_DF = pd.read_csv(ROOT / \"train.csv\") \nTEST_DF  = pd.read_csv(ROOT / \"sample_submission.csv\")\n\nTRAIN_DF.head(10)","metadata":{"execution":{"iopub.status.busy":"2022-08-16T12:38:36.688422Z","iopub.execute_input":"2022-08-16T12:38:36.689717Z","iopub.status.idle":"2022-08-16T12:38:36.784420Z","shell.execute_reply.started":"2022-08-16T12:38:36.689646Z","shell.execute_reply":"2022-08-16T12:38:36.783451Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def display_image(img, size=None, title=None):\n    pil_img = Img.fromarray(img.astype(np.uint8))\n    if size is not None:\n        pil_img = pil_img.resize((size, size))\n    display(pil_img)\n\ndef get_image(idx:int=None, name:str=None, split=\"train\", debug=True, path=ROOT):\n    df = TRAIN_DF if split=='train' else TEST_DF\n    label = -1\n    if idx == None:\n        idx = random.randrange(len(df))\n    if name is not None:\n        idx = int(df.loc[df[\"id\"] == name].index[0])\n    proposed_name, label = df.iloc[idx]\n    if name is None:\n        name = proposed_name\n    if split == 'test':\n        label = -1\n    filename = f\"{path}/{split}/{name}.jpeg\"\n    img = cv2.imread(filename, 0)\n    if debug:\n        print(idx, name, 'split', split, 'digit_sum', label)\n    return img, label\n\nimage, label = get_image()\ndisplay_image(image, 400)","metadata":{"execution":{"iopub.status.busy":"2022-08-16T12:38:36.785955Z","iopub.execute_input":"2022-08-16T12:38:36.786915Z","iopub.status.idle":"2022-08-16T12:38:37.020498Z","shell.execute_reply.started":"2022-08-16T12:38:36.786867Z","shell.execute_reply":"2022-08-16T12:38:37.019619Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Function to plot a list of angles","metadata":{}},{"cell_type":"code","source":"def draw_angles(angles, radius=1, deg=False):\n    if not isinstance(angles, (list, tuple)):\n        angles = list(angles)\n\n    # Draw circle\n    x0, y0 = 0, 0\n    circle = plt.Circle((x0, y0), radius, color='red', fill=False)\n    fig, ax = plt.subplots()\n    ax.add_artist(circle)\n    \n    for angle in angles:\n        if deg:\n            angle = angle * np.pi / 180\n        # Draw radius arrow\n        head_length = 0.05\n        ax.arrow(\n            x0, y0,\n            (radius - head_length) * np.cos(angle),\n            (radius - head_length) * np.sin(angle),\n            head_width=head_length, \n            head_length=head_length, \n            fc='k', ec='k'\n        )\n        \n    ax.set_xlim([-radius, radius])\n    ax.set_ylim([-radius, radius])\n    ax.set_aspect('equal')\n    ax.axis(False)\n    plt.show()\n    \ndraw_angles((45, 90), deg=True)","metadata":{"execution":{"iopub.status.busy":"2022-08-16T12:38:44.645053Z","iopub.execute_input":"2022-08-16T12:38:44.645353Z","iopub.status.idle":"2022-08-16T12:38:44.800980Z","shell.execute_reply.started":"2022-08-16T12:38:44.645321Z","shell.execute_reply":"2022-08-16T12:38:44.800011Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Detecting the grid !","metadata":{}},{"cell_type":"code","source":"def get_laplacian(img):\n    y = img\n    y = cv2.threshold(y, 100, 255, cv2.THRESH_BINARY)[1]\n    y = cv2.Laplacian(y, cv2.CV_64F)\n    y = cv2.convertScaleAbs(y)\n    return y\n\ndef process(img, rounded=0, debug=True):\n    RAD_2_DEG, DEG_2_RAD = 180 / np.pi, np.pi / 180\n    w, h = img.shape\n    x = img\n\n    angle_div = 180\n    # distance resolution in pixels of the Hough grid\n    # angular resolution in radians of the Hough grid\n    # 5  # minimum number of votes (intersections in Hough grid cell)\n    #1000 #50  # minimum number of pixels making up a line\n    # maximum gap in pixels between connectable line segments\n    rho = 1\n    theta = (np.pi / angle_div)\n    threshold = 50\n    min_line_length = 200 # 200 # 250 #500\n    max_line_gap = 50\n\n    \n    lines = cv2.HoughLinesP(x, rho, theta, threshold, np.array([]), min_line_length, max_line_gap)\n    if debug: print(lines.shape)\n    if lines is None:\n        return None\n    if debug: print('lines', lines.shape)\n    \n    # Run Hough on edge detected image\n    # Output \"lines\" is an array containing endpoints of detected line segments\n    out = np.zeros((w, h, 3), dtype=np.uint8)\n    neg = np.zeros((w, h), dtype=np.uint8)\n    A, D = [], []\n    for line in lines:\n        for x0,y0,x1,y1 in line:\n            dx,dy = x1-x0, y1-y0\n            if debug: assert dx>=0\n\n            d = np.sqrt((dx)**2+(dy)**2)\n            D.append(d)\n\n            a = np.rad2deg(-np.arctan2(dy, dx))\n            a = round(a, rounded)\n            A.append(a)\n             \n            cv2.line(out, (x0, y0), (x1, y1), (255, 0, 0), 5)\n            cv2.line(neg, (x0, y0), (x1, y1), 255, 15)\n                \n\n    if debug: display_image(out, 400)\n    A = np.array(A)\n    D = np.array(D)\n    \n    A_unique, A_counts = np.unique(A, return_counts=True)\n    if debug: print('DEBUG', A_unique)\n    if debug: print('DEBUG', A_counts)\n    E = A_counts / A_counts.sum()\n    #print('E', E)\n    W = np.array([D[a == A].sum() for a in A_unique]) / D.sum()\n    #print('W', W)\n    C = E * W\n    #print('C', C)\n    \n    # find top 2 angles\n    ind = np.argpartition(C, -2)[-2:]\n    if debug: print(ind, A_unique[ind])\n    if debug:draw_angles(A_unique, radius=1, deg=True)\n    if debug:draw_angles(A_unique[ind], radius=1, deg=True)\n    a_grid_0, a_grid_1 = A_unique[ind]\n    \n    \n    out = np.zeros((w, h, 3), dtype=np.uint8)\n    neg = np.zeros((w, h), dtype=np.uint8)\n    pos = np.zeros((w, h), dtype=np.uint8)\n    grid_lines = [[], []]\n    for line in lines:\n        for x0,y0,x1,y1 in line:\n            a = np.rad2deg(-np.arctan2(y1-y0, x1-x0))\n            a = round(a, rounded)\n            # GREEN - grid\n            \n            if any(np.abs(a-aref)<2 for aref in [a_grid_0]):\n                cv2.line(neg, (x0, y0), (x1, y1), 255, 10)\n                cv2.line(out, (x0, y0), (x1, y1), (0, 255, 0), 10)\n                grid_lines[0].append((x0,y0,x1,y1))\n                \n            elif any(np.abs(a-aref)<2 for aref in [a_grid_1]):\n                cv2.line(neg, (x0, y0), (x1, y1), 255, 10)\n                cv2.line(out, (x0, y0), (x1, y1), (0, 255, 0), 10)\n                grid_lines[1].append((x0,y0,x1,y1))\n                \n            elif any(np.abs(a-aref)<1 for aref in [0, 90, -90]): #exclude\n                cv2.line(neg, (x0, y0), (x1, y1), 255, 10)\n                cv2.line(out, (x0, y0), (x1, y1), (255, 0, 0), 10)\n            else:\n                cv2.line(pos, (x0, y0), (x1, y1),            255 ,  1)\n                cv2.line(out, (x0, y0), (x1, y1), (255, 255, 255), 10) # ignored\n\n    if debug: print('='*80)\n    if debug: display_image(out, 400)\n    #display_image(neg, 400)\n    \n    y = x.copy()\n    y[neg > 0] = 0\n    #y[pos > 0] = 255\n    if debug: display_image(y, 400)\n    return y, out, grid_lines\n\n\n#x, label = get_image(name='yrsmvijksq')\nx, label = get_image()\ndisplay_image(x, 400)\n\n\nprint('-'*80)\nx_lap = get_laplacian(x)\nx_r, grid, grid_lines = process(x_lap, debug=False)\ndisplay_image(x_r, 400)","metadata":{"execution":{"iopub.status.busy":"2022-08-16T12:51:45.512003Z","iopub.execute_input":"2022-08-16T12:51:45.512986Z","iopub.status.idle":"2022-08-16T12:51:46.410050Z","shell.execute_reply.started":"2022-08-16T12:51:45.512923Z","shell.execute_reply":"2022-08-16T12:51:46.408811Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"for _ in range(100):\n    x, label = get_image()\n    display_image(x, 400)\n\n\n    x_lap = get_laplacian(x)\n    x_r, grid, grid_lines = process(x_lap, debug=False)\n    display_image(x_r, 400)\n    \n    print('-'*80)\n","metadata":{"execution":{"iopub.status.busy":"2022-08-16T12:52:18.921390Z","iopub.execute_input":"2022-08-16T12:52:18.922330Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def det(a, b):\n        return a[0] * b[1] - a[1] * b[0]\n    \ndef line_intersection(line1, line2):\n    xdiff = (line1[0][0] - line1[1][0], line2[0][0] - line2[1][0])\n    ydiff = (line1[0][1] - line1[1][1], line2[0][1] - line2[1][1])\n    div = det(xdiff, ydiff)\n    if div == 0:\n        return None\n    d = (det(*line1), det(*line2))\n    x = det(d, xdiff) / div\n    y = det(d, ydiff) / div\n    return x, y\n\nS = 4000\np0, p1, p2, p3 = (0, 0), (0, S-1), (S-1, 0), (S-1, S-1)\nborders = [\n    (p0, p1), # left\n    (p1, p3), # bot\n    (p0, p2), # top\n    (p2, p3), # right\n]\n\ncolors = [(255, 0, 0), (0, 0, 255)]\n\n#y2 = np.zeros((S, S, 3), dtype=np.uint8)\ny2 = np.zeros((S, S), dtype=np.uint8)\ny = x.copy()\ny = cv2.cvtColor(y, cv2.COLOR_GRAY2RGB)\nfor i in range(2):\n    for (x0,y0,x1,y1) in grid_lines[i]:\n        a = (y1-y0)/(x1-x0)\n        b = y1-a*x1\n        line = ((x0,y0),(x1,y1))\n        \n        inters = [line_intersection(line, b) for b in borders]\n        inters = [tuple(map(round, p)) for p in inters if p is not None]\n        inters = [(x, y) for (x, y) in inters if (0<=x<=S) and (0<=y<=S)]\n        #print(inters)\n        cv2.line(y, inters[0], inters[1], colors[i], 5)\n        cv2.line(y2, inters[0], inters[1], random.randrange(255), 5)\n\n\ndisplay_image(y, 400)\ndisplay_image(y2, 400)","metadata":{"execution":{"iopub.status.busy":"2022-04-16T15:42:19.319084Z","iopub.execute_input":"2022-04-16T15:42:19.319436Z","iopub.status.idle":"2022-04-16T15:42:19.58122Z","shell.execute_reply.started":"2022-04-16T15:42:19.319414Z","shell.execute_reply":"2022-04-16T15:42:19.580248Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# TODO : ","metadata":{}},{"cell_type":"markdown","source":"process the last image to extract the grid parameters (2 not parallel lines)\n\nline : ($\\theta_0$, offset, spacing)\n\n![grid](https://i.imgur.com/VaeZyG7.png)","metadata":{}},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]}]}