Hough_with_DBSCAN_Preprocessing.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {
    "_cell_guid": "1bc74aa1-b9e7-43df-a128-2fd053a2f57a",
    "_uuid": "af970460c50b015f8d624d2566f0fe2545f0cf37"
   },
   "source": [
    "## Hough kernel with xyz-transform from DBSCAN kernel"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "_cell_guid": "16647aab-9afa-4213-a696-894f62b44de2",
    "_uuid": "f33a9129cc0f5d8366115ee14e857bb6819fe87b",
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "%matplotlib inline\n",
    "import matplotlib.pyplot as plt\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "import os\n",
    "from trackml.dataset import load_event, load_dataset\n",
    "from trackml.score import score_event\n",
    "from sklearn.preprocessing import StandardScaler\n",
    "\n",
    "path2data = '/home/alepfu/Desktop/dataset/train_100_events/'\n",
    "event_id = 'event000001000'\n",
    "hits, cells, particles, truth = load_event(path2data + event_id)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "_cell_guid": "eb7dd629-8b33-4ffa-9827-2edde2e2cb48",
    "_uuid": "9055e4033de48ec84f584cef1c95c9e78a704b7c",
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def cartesian_to_cylindrical(x, y, z):\n",
    "    \n",
    "    r = np.sqrt(x**2 + y**2)\n",
    "    phi = np.arctan2(y, x)\n",
    "    z = z\n",
    "    \n",
    "    return r, phi, z\n",
    "\n",
    "\n",
    "def create_hough_matrix(hits):\n",
    "    \n",
    "    # Preprocessing from DBSCAN kernel\n",
    "    hits_prep = hits[['hit_id', 'x', 'y', 'z']].copy()\n",
    "    hits_prep['r_xy'] = hits_prep.apply(lambda h: np.sqrt(h['x']**2 + h['y']**2), axis=1)\n",
    "    hits_prep['r_xyz'] = hits_prep.apply(lambda h: np.sqrt(h['x']**2 + h['y']**2 + h['z']**2), axis=1)\n",
    "    hits_prep['x2'] = hits_prep.apply(lambda h: h['x']/h['r_xyz'], axis=1)\n",
    "    hits_prep['y2'] = hits_prep.apply(lambda h: h['y']/h['r_xyz'], axis=1)\n",
    "    hits_prep['z2'] = hits_prep.apply(lambda h: h['z']/h['r_xy'], axis=1)\n",
    "    scaler = StandardScaler()\n",
    "    hits_prep[['x2', 'y2', 'z2']] = scaler.fit_transform(hits_prep[['x2', 'y2', 'z2']])\n",
    "    \n",
    "    hough_matrix = hits_prep[['hit_id', 'x2', 'y2', 'z2']].copy()\n",
    "    \n",
    "    hough_matrix.columns = ['HitID', 'X', 'Y', 'Z']\n",
    "    \n",
    "    r, phi, z = cartesian_to_cylindrical(hough_matrix['X'].values, hough_matrix['Y'].values, hough_matrix['Z'].values)\n",
    "\n",
    "    hough_matrix['R'] = r\n",
    "    hough_matrix['Phi'] = phi\n",
    "    \n",
    "    return hough_matrix\n",
    "\n",
    "\n",
    "def add_theta(hough_matrix, theta):\n",
    "    hough_matrix['Theta'] = theta\n",
    "    return hough_matrix\n",
    "\n",
    "def add_r0_inv(hough_matrix):\n",
    "    hough_matrix['R0Inv'] = (2. * np.cos(hough_matrix['Phi'] - hough_matrix['Theta']) / hough_matrix['R']).values\n",
    "    return hough_matrix\n",
    "\n",
    "def add_gamma(hough_matrix):\n",
    "    hough_matrix['Gamma'] = hough_matrix['Z']/hough_matrix['R']\n",
    "    return hough_matrix\n",
    "\n",
    "\n",
    "def digitize_column(hough_matrix, col, N, min_val=None, max_val=None):\n",
    "    \n",
    "    x = hough_matrix[col].values\n",
    "    if min_val is not None and max_val is not None:\n",
    "        bins = np.linspace(min_val, max_val, N)\n",
    "    else:\n",
    "        bins = np.linspace(x.min(), x.max(), N)\n",
    "    bin_ids = np.digitize(x, bins)\n",
    "    hough_matrix[col+'Digi'] = bin_ids\n",
    "    \n",
    "    return hough_matrix\n",
    "\n",
    "\n",
    "def combine_digi(hough_matrix, columns):\n",
    "    \n",
    "    hough_matrix['ComboDigi'] = np.zeros(len(hough_matrix))\n",
    "    \n",
    "    for i_col, acol in enumerate(columns):\n",
    "        digi = hough_matrix[acol]\n",
    "        hough_matrix['ComboDigi'] += digi * 10**(i_col * 5)\n",
    "    \n",
    "    return hough_matrix\n",
    "\n",
    "\n",
    "def count_combo_digi(hough_matrix):\n",
    "    \n",
    "    unique, indeces, counts = np.unique(hough_matrix['ComboDigi'].values, \n",
    "                                     return_counts=True, return_inverse=True)\n",
    "    hough_matrix['ComboDigiCounts'] = counts[indeces]\n",
    "    \n",
    "    return hough_matrix\n",
    "\n",
    "def out_of_border_counters_to_zero(hough_matrix, col, N):\n",
    "    hough_matrix['ComboDigiCounts'] *= (hough_matrix[col].values != 0) * (hough_matrix[col].values != N)\n",
    "    return hough_matrix\n",
    "\n",
    "def one_slice(hough_matrix, theta, N_bins_r0inv, N_bins_gamma, min_hits):\n",
    "        \n",
    "    tracks = []\n",
    "    \n",
    "    hough_matrix = add_theta(hough_matrix, theta)\n",
    "    hough_matrix = add_r0_inv(hough_matrix)\n",
    "    hough_matrix = add_gamma(hough_matrix)\n",
    "\n",
    "    hough_matrix = digitize_column(hough_matrix, 'R0Inv', N_bins_r0inv, None, None) # Set to None, to use min/max\n",
    "    hough_matrix = digitize_column(hough_matrix, 'Gamma', N_bins_gamma, None, None) \n",
    "\n",
    "    hough_matrix = combine_digi(hough_matrix, ['R0InvDigi', 'GammaDigi'])\n",
    "    hough_matrix = count_combo_digi(hough_matrix)\n",
    "\n",
    "    hough_matrix = out_of_border_counters_to_zero(hough_matrix, 'R0InvDigi', N_bins_r0inv)\n",
    "    hough_matrix = out_of_border_counters_to_zero(hough_matrix, 'GammaDigi', N_bins_gamma)\n",
    "    \n",
    "    counts = hough_matrix.ComboDigiCounts.values\n",
    "    bins = hough_matrix.ComboDigi.values\n",
    "    hit_ids = np.arange(len(hough_matrix))\n",
    "    for abin in np.unique(bins[counts >= min_hits]):\n",
    "        atrack = hit_ids[(bins == abin)]\n",
    "        tracks.append(atrack)\n",
    "        \n",
    "    return tracks, hough_matrix\n",
    "\n",
    "\n",
    "\n",
    "\n",
    "class Clusterer(object):\n",
    "    \n",
    "    def __init__(self, N_bins_r0inv, N_bins_gamma, N_theta, min_hits):\n",
    "        \n",
    "        self.N_bins_r0inv = N_bins_r0inv \n",
    "        self.N_bins_gamma = N_bins_gamma\n",
    "        self.N_theta = N_theta\n",
    "        self.min_hits = min_hits\n",
    "    \n",
    "    def predict(self, hits):\n",
    "        \n",
    "        tracks = []\n",
    "\n",
    "        hough_matrix = create_hough_matrix(hits)\n",
    "        for theta in np.linspace(-np.pi, np.pi, self.N_theta):\n",
    "            slice_tracks, hough_matrix = one_slice(hough_matrix, theta, self.N_bins_r0inv, self.N_bins_gamma, self.min_hits)\n",
    "            tracks += list(slice_tracks)\n",
    "\n",
    "        labels = np.zeros(len(hits))\n",
    "        used = np.zeros(len(hits))\n",
    "        track_id = 0\n",
    "        for atrack in tracks:\n",
    "            u_track = atrack[used[atrack] == 0]\n",
    "            if len(u_track) >= self.min_hits:\n",
    "                labels[u_track] = track_id\n",
    "                used[u_track] = 1\n",
    "                track_id += 1\n",
    "            \n",
    "        return labels"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "_cell_guid": "8c0c0c64-a5bc-4f4a-9e8c-ec29762c3627",
    "_uuid": "bea0436b431edad5b33e9974a74a82915eb0b2a0",
    "collapsed": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "CPU times: user 1min 51s, sys: 1.48 s, total: 1min 52s\n",
      "Wall time: 1min 50s\n"
     ]
    }
   ],
   "source": [
    "%%time\n",
    "# Warning: it takes about 100s per one event.\n",
    "\n",
    "model = Clusterer(N_bins_r0inv=200, N_bins_gamma=500, N_theta=500, min_hits=9)\n",
    "labels = model.predict(hits)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "_cell_guid": "b965aa18-af14-4e00-bae2-256b12b14276",
    "_uuid": "4f805ff897b36fac375740a73a877361ce69668f",
    "collapsed": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Your score:  0.0\n"
     ]
    }
   ],
   "source": [
    "def create_one_event_submission(event_id, hits, labels):\n",
    "    sub_data = np.column_stack(([event_id]*len(hits), hits.hit_id.values, labels))\n",
    "    submission = pd.DataFrame(data=sub_data, columns=[\"event_id\", \"hit_id\", \"track_id\"]).astype(int)\n",
    "    return submission\n",
    "\n",
    "submission = create_one_event_submission(0, hits, labels)\n",
    "score = score_event(truth, submission)\n",
    "\n",
    "print(\"Your score: \", score)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "_cell_guid": "e904efa3-7fc1-42ea-98b3-82d5011c93d8",
    "_uuid": "76f025c469e5964a6388dc138fafb3f875afdfad",
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.4"
  },
  "toc": {
   "nav_menu": {},
   "number_sections": false,
   "sideBar": true,
   "skip_h1_title": false,
   "toc_cell": false,
   "toc_position": {},
   "toc_section_display": "block",
   "toc_window_display": false
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}