{"cells":[{"metadata":{"_cell_guid":"d359b6c8-1803-4a89-888c-5273c501ead3","_uuid":"39eee1c9f51be754aa8282953e891b33f197b593"},"cell_type":"markdown","source":"![](https://upload.wikimedia.org/wikipedia/commons/1/1a/Tyre_Marks_in_the_Sand_-_geograph.org.uk_-_1546451.jpg)\n\n# Quick Problem Description\nLink every *track* to one *hit*.\n\n*Disclaimer: I'm oversimplifying physics to explain this problem*.\n\nEvery particle leaves a track behind it, like a car leaving tire marks in the sand.  We did not catch the particle in action.  Now we want to link every track (tire mark) to one hit that the particle created.\n\nIn every **event**, a large number of **particles** are released.  They move along a path leaving behind their **tracks**.  They eventually **hit** a particle detector surface on the other end.\n\nIn the training data we have the following information on each **event**:\n- **Hits**: $x, y, z$ coordinates of each hit on the particle detector\n- **Particles**: Each particle's initial position ($v_x, v_y, v_z$), momentum ($p_x, p_y, p_z$), charge ($q$) and number of hits\n- **Truth**: Mapping between hits and generating particles; the particle's trajectory, momentum and the hit weight\n- **Cells**: Precise location of where each particle hit the detector and how much energy it deposited\n"},{"metadata":{"_cell_guid":"089ad293-ca5d-4352-9e87-35300a757f3c","_uuid":"3d4c96fa4e9f4a55f70ac8fcf0e9073c18276757"},"cell_type":"markdown","source":"# Data Exploration:\n\n#### Import `trackml-library`\nThe easiest and best way to load the data is with the [trackml-library] that was built for this purpose.\n\nUnder your kernel's *Settings* tab -> *Add a custom package* -> *GitHub user/repo* (LAL/trackml-library)\n\nRestart your  kernel'\ns session and you will be good to go.\n\n[trackml-library]: https://github.com/LAL/trackml-library"},{"metadata":{"scrolled":true,"collapsed":true,"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","trusted":false},"cell_type":"code","source":"import os\n\nimport numpy as np\nimport pandas as pd\n\nfrom trackml.dataset import load_event\nfrom trackml.randomize import shuffle_hits\nfrom trackml.score import score_event\n\nimport matplotlib.pyplot as plt\nfrom mpl_toolkits.mplot3d import Axes3D\nimport seaborn as sns\n%matplotlib inline","execution_count":null,"outputs":[]},{"metadata":{"collapsed":true,"_uuid":"59c07c9273966e79b26594a4804f209064a94048","_cell_guid":"03e8af96-85ba-4d79-8601-453a32b8942e","trusted":false},"cell_type":"code","source":"event_prefix = 'event000001000'\nhits, cells, particles, truth = load_event(os.path.join('../input/train_1', event_prefix))\n\nmem_bytes = (hits.memory_usage(index=True).sum() \n             + cells.memory_usage(index=True).sum() \n             + particles.memory_usage(index=True).sum() \n             + truth.memory_usage(index=True).sum())\nprint('{} memory usage {:.2f} MB'.format(event_prefix, mem_bytes / 2**20))","execution_count":null,"outputs":[]},{"metadata":{"_cell_guid":"ad9c44f5-4b12-4aac-a557-d16fbbaa28be","_uuid":"967a443483052831a49bc0a019b7766c972a287e"},"cell_type":"markdown","source":"## Hits Data\n\n### Where Did it Hit?\nHere we have the $x, y, z$ global coordinates (in millimeters) of where the particles hit the detector surface."},{"metadata":{"_kg_hide-input":true,"collapsed":true,"_uuid":"bdaafeba9a5f961c74a486447729c0e24ae1c383","_cell_guid":"6a8fc52e-43c2-4de3-908c-506731b60e69","trusted":false},"cell_type":"code","source":"hits.head()","execution_count":null,"outputs":[]},{"metadata":{"_cell_guid":"1ab2bcf6-c621-4e20-a307-3e4f799271e2","_uuid":"f5cc5dff4d8c52263b0c11bf43a37f6e7177c053"},"cell_type":"markdown","source":"Here is the distribution of $x, y, z$ location of hits in event 1000.  This is only for one out of 8,850 events.\n\n### Vertical Intersection ($x, y$) in Detection Layers\nAs shown in the figure below, the hits are semi evenly distributed on the detector surface $x, y$.  The white circle in the center of the plot is where the beam pipe lies.  Thanks [agerom] for [the clarification][clar].\n\nThe colors represent different detector volumes.  Thanks to [Joshua Bonatt's notebook][josh].\n\n[josh]: https://www.kaggle.com/jbonatt/trackml-eda-etc\n[clar]: https://www.kaggle.com/wesamelshamy/trackml-problem-explanation-and-data-exploration/comments#323803\n[agerom]: https://www.kaggle.com/artemiosgeromitsos"},{"metadata":{"_kg_hide-input":true,"collapsed":true,"_uuid":"3fa89a0e36b380d73da92add20d844b6c44c2394","_cell_guid":"b9103439-c5b9-4f00-bc6f-c30186b37199","scrolled":false,"trusted":false},"cell_type":"code","source":"g = sns.jointplot(hits.x, hits.y,  s=1, size=12)\ng.ax_joint.cla()\nplt.sca(g.ax_joint)\n\nvolumes = hits.volume_id.unique()\nfor volume in volumes:\n    v = hits[hits.volume_id == volume]\n    plt.scatter(v.x, v.y, s=3, label='volume {}'.format(volume))\n\nplt.xlabel('X (mm)')\nplt.ylabel('Y (mm)')\nplt.legend()\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{"_cell_guid":"f5e39b13-889c-4b09-b3cd-49ff28d6a0db","_uuid":"8dce3b23a9c1e546d9d4a533bb5fa26207bde12e"},"cell_type":"markdown","source":"### Horizonal Intersection ($y, z$) in Detection Layers\nYou can think of the chart below as a horizontal intersection in the detection surface, where every dot is a hit.  Notice the relationship between the different activity levels in this chart and the one above for $x, y$.\n\nAgain, the colors represent different volumes in the detector surface."},{"metadata":{"_kg_hide-input":true,"collapsed":true,"_uuid":"01e21efd1d0fbe9ef58d6d0b108778b497e74418","_cell_guid":"2b8ebb84-6ba3-460b-8f0b-1cbc4c13bde8","trusted":false},"cell_type":"code","source":"g = sns.jointplot(hits.z, hits.y, s=1, size=12)\ng.ax_joint.cla()\nplt.sca(g.ax_joint)\n\nvolumes = hits.volume_id.unique()\nfor volume in volumes:\n    v = hits[hits.volume_id == volume]\n    plt.scatter(v.z, v.y, s=3, label='volume {}'.format(volume))\n\nplt.xlabel('Z (mm)')\nplt.ylabel('Y (mm)')\nplt.legend()\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{"_cell_guid":"8317678a-0b46-4aae-ab87-63f2a8f34627","_uuid":"e810d9f26e7098ee93d2a759c0e5cd3a3d871b37"},"cell_type":"markdown","source":"And here is how the hits in this event look like in 3D.  Again, a sample from one event.  This combines the previous two charts in 3D.\n\nNotice how the particles penetrate the detector surface along $z$ coordinate."},{"metadata":{"_kg_hide-input":true,"collapsed":true,"_uuid":"8b12e47ba2b0390593de00b5831f2967273d2727","_cell_guid":"5a459a75-7177-486e-9672-5d83c57c9b32","trusted":false},"cell_type":"code","source":"fig = plt.figure(figsize=(12, 12))\nax = fig.add_subplot(111, projection='3d')\nfor volume in volumes:\n    v = hits[hits.volume_id == volume]\n    ax.scatter(v.z, v.x, v.y, s=1, label='volume {}'.format(volume), alpha=0.5)\nax.set_title('Hit Locations')\nax.set_xlabel('Z (millimeters)')\nax.set_ylabel('X (millimeters)')\nax.set_zlabel('Y (millimeters)')\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{"_cell_guid":"d900418a-d5a0-47de-b38f-4d951674436f","_uuid":"85000e23ede450d2ed202ac2e90885a51e59a7e6"},"cell_type":"markdown","source":"### Affected Surface Object\nThe **volume**, **layer** and **module** are nested parts on the detector surface.  The volume is made of layers, which in turn have modules.  Analyzing their response could help us understand if some of them are dead/defective and therefore we may need to account for the bias they cause.\n\nThe figure below shows a plot of every combination of `x`, `y`, `volume`, `layer` and `module`.  The colors identify different *volumes*.  Along the main diagonal we have the variables' histograms.\n\nThe (`hit_id`, `x`) and (`hit_id`, `y`) pairs show us how different volumes are layered.\n\nThis figure idea is taken from [Joshua Bonatt's notebook][josh].\n\n[josh]: https://www.kaggle.com/jbonatt/trackml-eda-etc"},{"metadata":{"_kg_hide-input":true,"collapsed":true,"_uuid":"220569de9ac9b768ed21eb45fce93f29cf715301","_cell_guid":"57008c54-dd96-4078-b15c-4309df70fcf1","scrolled":false,"trusted":false},"cell_type":"code","source":"hits_sample = hits.sample(8000)\nsns.pairplot(hits_sample, hue='volume_id', size=8)\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{"_cell_guid":"9d1f9269-5cba-47f7-8c63-232217e14dfa","_uuid":"4cb917ec90ef3557fd74aefec7add41f229e541c"},"cell_type":"markdown","source":"## Particle Data\nThe particle data help us understand each particle's initial position, momentum, and charge, which we can join with the event truth data set to get the particle's final position and momentum.  This is needed to identify the tracks that each particle generated.\n\nThe data look like this:"},{"metadata":{"_kg_hide-input":true,"collapsed":true,"_uuid":"fa553f014a7e9cd86fdd63ca6f4cf3e509ed77ee","_cell_guid":"767168f3-72a8-4a18-acdd-d2f99cec82ca","scrolled":true,"trusted":false},"cell_type":"code","source":"particles.head()","execution_count":null,"outputs":[]},{"metadata":{"_cell_guid":"b4609d74-9f85-4864-853a-5d7ba6648a86","_uuid":"3b158a527abd9fb972b2c93e61eb723b9584cd7b"},"cell_type":"markdown","source":"### Hit Rate and Charge Distribution\nLet's see the distribution of the number of hits per particle, shown below.  A significant number of particles had no attributed hits, and most of them have positive charge in this event."},{"metadata":{"_kg_hide-input":true,"collapsed":true,"_uuid":"6d01bb62a80675ab2bcf6d7b4e77138cdbef2c68","_cell_guid":"8e3fe657-91af-41ac-9103-c1fbab7b87bc","scrolled":false,"trusted":false},"cell_type":"code","source":"plt.figure(figsize=(15, 5))\nplt.subplot(1, 2, 1)\nsns.distplot(particles.nhits.values, axlabel='Hits/Particle', bins=50)\nplt.title('Distribution of number of hits per particle for event 1000.')\nplt.subplot(1, 2, 2)\nplt.pie(particles.groupby('q')['vx'].count(),\n        labels=['negative', 'positive'],\n        autopct='%.0f%%',\n        shadow=True,\n        radius=0.8)\nplt.title('Distribution of particle charges.')\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{"_cell_guid":"2758c4e7-6f42-4231-bae5-eed457b5a0a9","_uuid":"59debf3c9e7b418f2fda294c17efeb177d5bb554"},"cell_type":"markdown","source":"### Initial Position and Momentum\nLet's now take a look at the initial position of the particles around the global coordinates' origin $(x, y)=(0,0)$, as shown in the figure below.\n\nThe initial position distribution is more concentrated around the origin (less variance) than its hit position (shown above under the Hits Data section).  As the particles hit the detection surface, they tend to scatter as shown in the particle trajectory plot at the end of this notebook.\n\nThe colors here show the number of hits for each particle."},{"metadata":{"_kg_hide-input":true,"collapsed":true,"_uuid":"00bfca1a73f34b3e57c07adb16085f4dac54d4d4","_cell_guid":"ac738872-a844-4689-80fd-7a890f4ceee6","trusted":false},"cell_type":"code","source":"g = sns.jointplot(particles.vx, particles.vy,  s=3, size=12)\ng.ax_joint.cla()\nplt.sca(g.ax_joint)\n\nn_hits = particles.nhits.unique()\nfor n_hit in n_hits:\n    p = particles[particles.nhits == n_hit]\n    plt.scatter(p.vx, p.vy, s=3, label='Hits {}'.format(n_hit))\n\nplt.xlabel('X (mm)')\nplt.ylabel('Y (mm)')\nplt.legend()\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{"_cell_guid":"e2f2ef2c-29c2-4f74-bbd3-f013bc9feebd","_uuid":"42c65ac6fc47295acb8dd53263ebdd5d13bc0b59"},"cell_type":"markdown","source":"And here is the initial position of the particles in a $z$, $y$ view.  Colors show number of hits."},{"metadata":{"_kg_hide-input":true,"collapsed":true,"_uuid":"7d9f003a21cd7a6a90aa7c836902861d8403f2fa","_cell_guid":"efd8b505-251d-46a4-9dd8-ceda7b681cf3","trusted":false},"cell_type":"code","source":"g = sns.jointplot(particles.vz, particles.vy,  s=3, size=12)\ng.ax_joint.cla()\nplt.sca(g.ax_joint)\n\nn_hits = particles.nhits.unique()\nfor n_hit in n_hits:\n    p = particles[particles.nhits == n_hit]\n    plt.scatter(p.vz, p.vy, s=3, label='Hits {}'.format(n_hit))\n\nplt.xlabel('Z (mm)')\nplt.ylabel('Y (mm)')\nplt.legend()\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{"_cell_guid":"59445b5a-dbd6-4172-861c-9c030697eef5","_uuid":"9d9985b942181ca712f7868e1ff8526677f2717d"},"cell_type":"markdown","source":"And this is what they look like in 3D."},{"metadata":{"_kg_hide-input":true,"collapsed":true,"_uuid":"ebf24188198254c2b7b2061596d6292430c26233","_cell_guid":"89143b53-2a3e-430a-969e-b9573fd1f65e","trusted":false},"cell_type":"code","source":"fig = plt.figure(figsize=(12, 12))\nax = fig.add_subplot(111, projection='3d')\nfor charge in [-1, 1]:\n    q = particles[particles.q == charge]\n    ax.scatter(q.vz, q.vx, q.vy, s=1, label='Charge {}'.format(charge), alpha=0.5)\nax.set_title('Sample of 1000 Particle initial location')\nax.set_xlabel('Z (millimeters)')\nax.set_ylabel('X (millimeters)')\nax.set_zlabel('Y (millimeters)')\nax.legend()\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{"_cell_guid":"6cc06ab2-1953-4070-a65d-0b6222254e2c","_uuid":"b64e8f430428799b25c882457d38c74ab65fbf70"},"cell_type":"markdown","source":"## Pair plot\nLet's now take a look at the relationship between different pair combinations of the particle variables.  Again, the colors represent the number of hits.\n\nThere is no large skew in the distribution of the number of hits over other variables.  It looks like the particles are targetted towards the global origin $(x,y)=(0,0)$ and are evenly distributed around it."},{"metadata":{"_kg_hide-input":true,"collapsed":true,"_uuid":"275576a476f7a5f9fbcab8706fa40ca7280b6ccc","_cell_guid":"e94570e5-a40d-46cc-841f-177ccd8eef5b","trusted":false},"cell_type":"code","source":"p_sample = particles.sample(8000)\nsns.pairplot(p_sample, vars=['particle_id', 'vx', 'vy', 'vz', 'px', 'py', 'pz', 'nhits'], hue='nhits', size=8)\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{"_cell_guid":"cc766f4f-b28e-44b0-82f3-12159e6a6687","_uuid":"26cec5efdc80f7cf931b5b888e014ee83633edb9"},"cell_type":"markdown","source":"### Particle Trajectory\nWe can reconstruct the trajectories for a few particles given their intersection points with the detection layers.  As explained in the [competition evaluation page][ceval], hits from straight tracks have larger wieghts, and random tracks or hits with very short tracks have weights of zero.  The figure below shows two such examples.\n\nThanks to [maka's notebook][traj] for the idea.\n\n[traj]: https://www.kaggle.com/makahana/quick-trajectory-plot\n[ceval]: https://www.kaggle.com/c/trackml-particle-identification#evaluation"},{"metadata":{"_kg_hide-input":true,"collapsed":true,"_uuid":"f80be8052d98f82af6bc17a3726a792ff36f892c","_cell_guid":"bfd1e6ec-ed44-4ece-bd54-51f60c6cb7b8","trusted":false},"cell_type":"code","source":"# Get particle id with max number of hits in this event\nparticle = particles.loc[particles.nhits == particles.nhits.max()].iloc[0]\nparticle2 = particles.loc[particles.nhits == particles.nhits.max()].iloc[1]\n\n# Get points where the same particle intersected subsequent layers of the observation material\np_traj_surface = truth[truth.particle_id == particle.particle_id][['tx', 'ty', 'tz']]\np_traj_surface2 = truth[truth.particle_id == particle2.particle_id][['tx', 'ty', 'tz']]\n\np_traj = (p_traj_surface\n          .append({'tx': particle.vx, 'ty': particle.vy, 'tz': particle.vz}, ignore_index=True)\n          .sort_values(by='tz'))\np_traj2 = (p_traj_surface2\n          .append({'tx': particle2.vx, 'ty': particle2.vy, 'tz': particle2.vz}, ignore_index=True)\n          .sort_values(by='tz'))\n\nfig = plt.figure(figsize=(10, 10))\nax = fig.add_subplot(111, projection='3d')\n\nax.plot(\n    xs=p_traj.tx,\n    ys=p_traj.ty,\n    zs=p_traj.tz, marker='o')\nax.plot(\n    xs=p_traj2.tx,\n    ys=p_traj2.ty,\n    zs=p_traj2.tz, marker='o')\n\nax.set_xlabel('X (mm)')\nax.set_ylabel('Y (mm)')\nax.set_zlabel('Z  (mm) -- Detection layers')\nplt.title('Trajectories of two particles as they cross the detection surface ($Z$ axis).')\nplt.show()","execution_count":null,"outputs":[]}],"metadata":{"language_info":{"name":"python","version":"3.6.5","mimetype":"text/x-python","codemirror_mode":{"name":"ipython","version":3},"pygments_lexer":"ipython3","nbconvert_exporter":"python","file_extension":".py"},"kernelspec":{"display_name":"Python 3","language":"python","name":"python3"}},"nbformat":4,"nbformat_minor":1}