{"cells":[{"metadata":{"collapsed":true,"_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5"},"cell_type":"markdown","source":"# TrackML Particle Tracking Challenge"},{"metadata":{"collapsed":true,"_cell_guid":"79c7e3d0-c299-4dcb-8224-4455121ee9b0","_uuid":"d629ff2d2480ee46fbb7e2d37f6b5fab8052498a"},"cell_type":"markdown","source":"This kernel is a code along of [Joshua Bonatt's TrackML EDA, etc.](https://www.kaggle.com/jbonatt/trackml-eda-etc/notebook) kernel. I want to get comfortable with reproducing his results more/less from scratch."},{"metadata":{"_cell_guid":"6033f185-b296-4a05-b382-3a9cdcb52031","_uuid":"f5220d5ba1ccb08eb1f25594a5f9d99e426c65de"},"cell_type":"markdown","source":"\n# Contents\n1. [Imports & Setup](#imports)\n2. [Hits](#hits)\n3. [Cells](#cells)\n4. [Particles](#particles)\n5. [Truth](#truth)"},{"metadata":{"_cell_guid":"a57cd087-1ffa-4378-b817-ff0cb9cbc75a","_uuid":"fd9536b2f5d67376ac8355b9dddfe601a6367e13"},"cell_type":"markdown","source":"## <a name=\"imports\">Imports & Setup</a>"},{"metadata":{"_cell_guid":"f59f1163-f4ad-475b-9ca9-6ed3894ad215","_uuid":"8b122b4d66155745f45cda4e57b2ee9455e6de2a"},"cell_type":"markdown","source":"[TrackML library](https://github.com/LAL/trackml-library) downloaded via: **Settings** >> **Add a custom package** >> *Github user/repo:* `LAL/trackml-library`"},{"metadata":{"collapsed":true,"_cell_guid":"bf89d301-f22f-4d66-9b07-aab482f9d1c0","_uuid":"313a06d790d500c78bf7e3a4cb3741a3067cdbb1","trusted":true},"cell_type":"code","source":"%matplotlib inline\n\nimport numpy as np\nimport pandas as pd\nimport matplotlib.pyplot as plt\n\nfrom mpl_toolkits import mplot3d\nimport seaborn as sns","execution_count":1,"outputs":[]},{"metadata":{"collapsed":true,"_cell_guid":"70432189-2103-4fd8-9e05-d739553c4bf1","_uuid":"1f21dc0508d86638d094f47077c05c7a69e6130c","trusted":true},"cell_type":"code","source":"from trackml.dataset import load_event, load_dataset\nfrom trackml.randomize import shuffle_hits\nfrom trackml.score import score_event","execution_count":2,"outputs":[]},{"metadata":{"_cell_guid":"12b48956-0ab8-4f03-afe8-d62d7bb7ceed","_uuid":"0ceadf7dca7dc924a96d013835dae56946e23659"},"cell_type":"markdown","source":"## Load one event for EDA"},{"metadata":{"collapsed":true,"_cell_guid":"56b35bb9-6e1a-4088-b85b-5be063f1d369","_uuid":"2f8fc090e06805a2f94a3bd1f9177ebed59dec4a","trusted":true},"cell_type":"code","source":"# One event of 8850\nevent_id = 'event000001000'\n# \"All methods either take or return pandas.DataFrame objects\"\nhits, cells, particles, truth = load_event('../input/train_1/' + event_id)","execution_count":3,"outputs":[]},{"metadata":{"_cell_guid":"e90da0f0-82fa-470d-b286-dc239dc8800b","_uuid":"33617f5ddc23d5551cc9e09a2a345c91dceda9c7"},"cell_type":"markdown","source":"## Data Files"},{"metadata":{"_cell_guid":"e7de255d-df27-4851-94be-aec8b0b67527","_uuid":"d7448e1a2698f94b294c44d090fae3a34935c4a0"},"cell_type":"markdown","source":"**hits**:\n- `hit_id`: numerical identifier of the hit inside the event\n- `x,y,z`: measured x,y,z position [mm] of hit in global coords\n- `volume_id`: numerical identifier of detector group\n- `layer_id`: numerical identifier of detector layer inside group\n- `module_id`: numerical identifer of detector module inside layer\n\n**cells**:\n- `hit_id`: num id of hit as defined in hits file\n- `ch0, ch1`: channel id/coords unique w/n 1 module\n- `value`: signal value information: particle charge deposition\n\n**particles**:\n- `particle_id`: num id of particle inside event\n- `vx,vy,vz`: inittial position or vertex [mm] in global coords\n- `px,py,pz`: intitial momentum [GeV/c] along each global axis\n- `q`: particle charge (multiple of absolute electron charge)\n- `nhits`: number of hits by this particle\n\n**truth**:\n- `hit_id`: num id of hit as defined in hits file\n- `particle_id`: num id of particle as defnd in particles file.\n    - `0`: hit did *not* originate from a reconstrible particle, ie: detector noise\n- `tx,ty,tz`: true intersection point in global coords [mm] between particle trajectory and sensing surface.\n- `tpx,tpy,tpz`: true particle momentum [GeV/c] in global coords at intersection point. Corresponding vector is *tangent* to particle traj at intersection point\n- `weight`: per-hit weight used for scoring metric. **`Σ`**`(weights in 1 event) = 1`\n"},{"metadata":{"_cell_guid":"188b1895-c695-4256-84e3-55373df4b921","_uuid":"08a0b72a6ed174eac5ef9d3b865f9b3296d641fd"},"cell_type":"markdown","source":"## <a name=\"hits\">Hits</a>"},{"metadata":{"_cell_guid":"4f1f1a07-c58a-48dd-8bf6-68b2146f24df","_uuid":"7d4d6b7e568a1fbbb5c56c8f3e06f1a85c75e273","trusted":true},"cell_type":"code","source":"hits.head()","execution_count":4,"outputs":[]},{"metadata":{"_cell_guid":"6747d972-de0a-46c6-b2ce-8cb82ef7d6b7","_uuid":"7faa649e20d22d3b2f1c95641c5f542f11d51f3b","trusted":true},"cell_type":"code","source":"hits.tail()","execution_count":5,"outputs":[]},{"metadata":{"_cell_guid":"e0ba0b5a-3783-42a1-b9e9-6e071c114179","_uuid":"998136bedc3bf071a933986e5ad04a41dacf440d","trusted":true},"cell_type":"code","source":"hits.describe()","execution_count":6,"outputs":[]},{"metadata":{"_cell_guid":"fffa674a-e383-4275-9f0a-e5d98dc865da","_uuid":"f8ca749f4e546ffc79089183c96e6452d018ce32"},"cell_type":"markdown","source":"The mean of `x,y,z` is only a few mm from the detector's center. The std is very high though: 305, 305, 1061. This means there's a lot of spread in hit location."},{"metadata":{"_cell_guid":"9ea2fdef-26c8-49fd-929e-90d051c2e30b","_uuid":"a900fcb3b9771bb48684efa16a1ea763b6df9162"},"cell_type":"markdown","source":"## Spatial Distribution of Hits"},{"metadata":{"_cell_guid":"80979124-3bdc-43c7-b958-84fc63972520","_uuid":"91610ec099b2f95e72538210ded26b84b53a0613","trusted":true},"cell_type":"code","source":"# plt.figure(figsize=(10,10))\n# plt.scatter(hits.x,hits.y, s=1)\n# plt.show()\n# Same as above, but include Univariate plots & Pearson correlation coeffs:\nradialview = sns.jointplot(hits.x, hits.y, size=10, s=1)\nradialview.set_axis_labels('x [mm]', 'y [mm]')\nplt.show()","execution_count":7,"outputs":[]},{"metadata":{"_cell_guid":"25cf09d2-7903-4cf9-8610-f20614f8daf3","_uuid":"f74d1fc42a22c915c38bb4ab262f60890c564dc9"},"cell_type":"markdown","source":"The solid core is likely more detectors. Zooming in to find out:"},{"metadata":{"_cell_guid":"315e680f-8de7-4a82-93a4-709b293e5a29","_uuid":"265548adebcc9471ebd6d71fb1b58e02ca02055d","trusted":true},"cell_type":"code","source":"radialview = sns.jointplot(hits[hits.x.abs() < 200].x, hits[hits.y.abs() < 200].y, size=10, s=1)\nradialview.set_axis_labels('x [mm]', 'y [mm]');\n# plt.show()","execution_count":8,"outputs":[]},{"metadata":{"_cell_guid":"c2fb24bb-be00-498d-8804-62c83a853254","_uuid":"b38b17b77d60e32b3ba3bf442277a1718f1a792c"},"cell_type":"markdown","source":"The scattering between rings are events from vertical detectors. Below shows these caps removed and again shows the concentric nature of the inner detector (*I still don't exactly know what that means -- I guess the scattering is further along (larger values) of hits on the z-detector*):\n\n$\\rightarrow$ *yeah, so you have cylindrical detectors near center, and discal detectors further out. Setting the z limit to 200 effectively ignores detections from the disks further away. So you only see discrete radial detections, instead of a spread.*"},{"metadata":{"collapsed":true,"_cell_guid":"1630771c-933c-425d-92fc-8f69a97896b4","_uuid":"e69058608d08a96f0922dc6fb4cb2fa9e3472c0b","trusted":true},"cell_type":"code","source":"def radial_display(dat, lim=2000):\n    radialview = sns.jointplot(dat[dat.x.abs()<lim].x, dat[dat.y.abs()<lim].y, size=10, s=1)\n    radialview.set_axis_labels('x [mm]', 'y [mm]')\n    plt.show()","execution_count":9,"outputs":[]},{"metadata":{"_cell_guid":"d5ed9979-b72a-4894-90d0-cba477cb3143","_uuid":"420b5429afb767a53c2d344910058f73d1a4913d","trusted":true},"cell_type":"code","source":"nocap = hits[hits.z.abs() < 200]\nradial_display(nocap, lim=200)","execution_count":10,"outputs":[]},{"metadata":{"_cell_guid":"59f6f65b-1df6-4280-9bc4-41305cf020a2","_uuid":"7e634335107e2ccd38b227fa9eec7107169521c1"},"cell_type":"markdown","source":"The detectors are layered as flat shingled rectangles. Zooming into the center-most detectors:"},{"metadata":{"_cell_guid":"e284f553-9b19-4a91-9216-343a8e9a3b6a","_uuid":"4440c0e4e34911e21475ff7af68f19cf5bab6003","trusted":true},"cell_type":"code","source":"radial_display(nocap, lim=50)","execution_count":11,"outputs":[]},{"metadata":{"_cell_guid":"4424a639-54da-4e11-bf88-42ba8589e34a","_uuid":"8a25681eaee23ae3e58b18f7a2192b9576a76936"},"cell_type":"markdown","source":""},{"metadata":{"_cell_guid":"3c8996b1-6f75-4bec-add5-02949708eb5d","_uuid":"b153e01ee262baaccddc2dda1b2ab1cfddda64e3","trusted":true},"cell_type":"code","source":"def side_display(fgsz = (24,8)): \n    plt.figure(figsize=fgsz)\n    axialview = plt.scatter(hits.z, hits.y, s=1)\n    plt.xlabel('z (mm)')\n    plt.ylabel('y (mm)')\n    plt.show()\nside_display()","execution_count":12,"outputs":[]},{"metadata":{"_cell_guid":"1b9ce5d5-6351-4ff5-9974-74d12b55fc78","_uuid":"286050233096ec4644f88fa8adcce54d017a4e42"},"cell_type":"markdown","source":"3D plot of a random sample of hits."},{"metadata":{"_cell_guid":"a95996b1-8ada-49ff-a920-850e9a830fd4","_uuid":"95ce098392f749d634e3bae1075a38ff5538e09c","trusted":true},"cell_type":"code","source":"def iso_display():\n    plt.figure(figsize=(15,15))\n    ax = plt.axes(projection='3d')\n    sample = hits.sample(30000)\n    ax.scatter(sample.z, sample.x, sample.y, s=5, alpha=0.5)\n    ax.set_xlabel('z (mm)')\n    ax.set_ylabel('x (mm)')\n    ax.set_zlabel('y (mm)')\n    # These two added to widen the 3D space\n    ax.scatter(3000,3000,3000, s=0)\n    ax.scatter(-3000,-3000,-3000, s=0)\n    plt.show()\niso_display()","execution_count":13,"outputs":[]},{"metadata":{"_cell_guid":"26a225fa-d006-4bbd-bee7-c3391cc0ad3e","_uuid":"bf0571faf4fc29a27fbb1a988c9814a80538472f"},"cell_type":"markdown","source":"## Location of Individual Detector Groups\n\nplotting each detector group as a different color:"},{"metadata":{"_cell_guid":"84e394b0-2eb9-4959-b022-ff7fff46d74f","_uuid":"80fde9d8a94c66e5855d79e721a8ac846a62bca0","trusted":true},"cell_type":"code","source":"volumes = hits.volume_id.unique()\n\nfg,ax = plt.subplots(figsize=(15,15))\nfor volume in volumes:\n    v = hits[hits.volume_id == volume]\n    ax.scatter(v.x, v.y, s=10, label='Volume '+str(volume), alpha=0.5)\nax.set_title('Detector Volumes, Radial View')\nax.set_xlabel('x [mm]')\nax.set_ylabel('y [mm]')\nax.legend()\nplt.show()\n\nvolumes = hits.volume_id.unique()\n\nfg,ax = plt.subplots(figsize=(24,8))\nfor volume in volumes:\n    v = hits[hits.volume_id == volume]\n    ax.scatter(v.z, v.y, s=10, label='Volume '+str(volume), alpha=0.5)\nax.set_title('Detector Volumes, Axial View')\nax.set_xlabel('z [mm]')\nax.set_ylabel('y [mm]')\nax.legend()\nplt.show()","execution_count":14,"outputs":[]},{"metadata":{"_cell_guid":"a921b132-2add-4004-8670-cc96c5a7184b","_uuid":"fc7eac6af38c3ff4b007ef6284aed7eca9c1505d"},"cell_type":"markdown","source":"Plotting this in 3D:"},{"metadata":{"_cell_guid":"9ff56632-40d1-4bf2-be17-2b82e5c62e88","_uuid":"94373a3202c9d7ba08f7347823bafc0c68b29ed0","trusted":true},"cell_type":"code","source":"sample = hits.sample(30000)\nplt.figure(figsize=(20,20))\nax = plt.axes(projection='3d')\nfor volume in volumes:\n    v = sample[sample.volume_id == volume]\n    ax.scatter(v.z, v.x, v.y, s=5, label='Volume '+str(volume), alpha=0.5)\nax.set_xlabel('z (mm)'); ax.set_ylabel('x (mm)'); ax.set_zlabel('y (mm)')\nax.legend()\n# added to widen the 3D space:\nax.scatter(3000,3000,3000, s=0); ax.scatter(-3000,-3000,-3000, s=0)\nplt.show()","execution_count":15,"outputs":[]},{"metadata":{"_cell_guid":"b744dc0e-f2ac-4975-96b1-0aa9aae5bdfe","_uuid":"d4f0df2b37aa94c2d890a4dcee399759ee917e5e"},"cell_type":"markdown","source":"We can also look at the layers:"},{"metadata":{"_cell_guid":"e7b42566-05e7-480d-b3ef-2f9bf035b7e1","_uuid":"7a26400fe7aad391c53de81ed89f6907e81d7e74","trusted":true,"scrolled":false},"cell_type":"code","source":"# RADIAL\nlayers = hits.layer_id.unique()\nfg,ax  = plt.subplots(figsize=(15,15))\nfor l_name in layers:\n    l = hits[hits.layer_id == l_name]\n    ax.scatter(l.x, l.y, s=10, label='Layer '+str(l_name), alpha=0.5)\nax.set_title('Detector Layers, Radial View'); ax.set_xlabel('x [mm]'); ax.set_ylabel('y [mm]')\nax.legend()\nplt.show()\n# AXIAL\nfg,ax  = plt.subplots(figsize=(24,8))\nfor l_name in layers:\n    l = hits[hits.layer_id == l_name]\n    ax.scatter(l.z, l.y, s=10, label='Layer '+str(l_name), alpha=0.5)\nax.set_title('Detector Layers, Axial View'); ax.set_xlabel('z [mm]'); ax.set_ylabel('y [mm]')\nax.legend()\nplt.show()\n# ISOMETRIC\nsample = hits.sample(30000)\nplt.figure(figsize=(20,20))\nax = plt.axes(projection='3d')\nfor layer in layers:\n    l = sample[sample.layer_id == layer]\n    ax.scatter(l.z, l.x, l.y, s=5, label='Layer '+str(layer), alpha=0.5)\nax.set_xlabel('z (mm)'); ax.set_ylabel('x (mm)'); ax.set_zlabel('y (mm)'); ax.legend()\n# added to widen the 3D space\nax.scatter(3000,3000,3000, s=0); ax.scatter(-3000,-3000,-3000, s=0)\nplt.show()","execution_count":21,"outputs":[]},{"metadata":{"collapsed":true,"_cell_guid":"afb051ef-c950-4114-bdbd-5d8467538c14","_uuid":"74943c1c9831e36f4bcfae1c795a5e890b4f6cba","trusted":false},"cell_type":"markdown","source":"modules not plotted for now (too many)"},{"metadata":{"collapsed":true,"_cell_guid":"e2f6df9e-d397-49d6-bc71-976cac2abf3b","_uuid":"33060f7ee1cf4484a4b4525b213b13342fe37cdb","trusted":false},"cell_type":"markdown","source":"## Detector Group Inquiry"},{"metadata":{"collapsed":true,"_cell_guid":"56c256e8-b1c1-40ce-b2f8-d76c823612c4","_uuid":"6688914b3ac5780fe99e90dd74bb126eed4962d8","trusted":false},"cell_type":"markdown","source":"Modules make up Layers. Layers make up Volumes. `module_id` is a subdir of `layer_id` which is a subdir of `volume_id`. Cells are the smallest unit of resolution and this a subdir of `module_id`. We can look at the population of these:"},{"metadata":{"_cell_guid":"5596dc61-0f73-49ba-ac3c-f54842237b0a","_uuid":"8c1d435547e304fd28b4bb3aca4d5ba37061be41","trusted":true},"cell_type":"code","source":"groups = [hits.volume_id, hits.layer_id, hits.module_id, cells.ch0, cells.ch1]\nfig,axes = plt.subplots(1,5, figsize=(30,10))\nfor i,ax in enumerate(axes):\n    sns.distplot(groups[i], ax=ax)","execution_count":38,"outputs":[]},{"metadata":{"collapsed":true,"_cell_guid":"dd2ccb0d-8381-49bc-a50a-bf3945a2431d","_uuid":"e0ab054a1a011c92103719499b1c8b09d5cea68c","trusted":false},"cell_type":"markdown","source":"A guess is that low-ID layers, modules, cells are closer to the center (due to the higher number of hits). We can plot hits by their radius to see clearer:"},{"metadata":{"_cell_guid":"56608e83-f817-4db4-ab06-2044b6462145","_uuid":"4030e509b9626842d67a605dd1c3fd907d2bfd1c","trusted":true},"cell_type":"code","source":"radius2 = np.sqrt(hits.x**2 + hits.y**2)\nradius3 = np.sqrt(hits.x**2 + hits.y**2 + hits.z**2)\nz2 = hits.z**2\nrads = [radius2, radius3, z2]\n\naxlbls = ['sqrt(x^2 + y^2)', 'sqrt(x^2 + y^2 + z^2)', 'z']\n\nfig,axes = plt.subplots(1,3, figsize=(30,10))\nfor i,ax in enumerate(axes):\n    sns.distplot(rads[i], axlabel=axlbls[i], ax=ax)","execution_count":41,"outputs":[]},{"metadata":{"collapsed":true,"_cell_guid":"5553403a-cadc-4525-a97b-ca297918a514","_uuid":"2be81cffbfa9841a15eeb9c1eb0e9ac2773127ce","trusted":false},"cell_type":"markdown","source":"The general distribution of events are proportional to the radius. Plotting groups by radius:"},{"metadata":{"_cell_guid":"4d00ddbf-fe92-407b-9e9d-17c1182ec924","_uuid":"899dd8716c2ec4c546a7dd3a082e4c132e01568a","trusted":true},"cell_type":"code","source":"labels = [['volume_id','radius'],['layer_id','radius'],['module_id','radius']]\ngroups = [hits.volume_id, hits.layer_id, hits.module_id]\nfig,axes = plt.subplots(1,3, figsize=(30,10))\nfor i,ax in enumerate(axes):\n    ax.set_xlabel(labels[i][0]); ax.set_ylabel(labels[i][1])\n    ax.scatter(groups[i], radius2)","execution_count":44,"outputs":[]},{"metadata":{"collapsed":true,"_cell_guid":"365c451a-6c80-440a-a24e-1ccdae31215c","_uuid":"cdf9d4851fc04fe95d00b77451155056abbeed9a","trusted":false},"cell_type":"markdown","source":"From these plots:\n- **Volumes** are named for **left**, **center**, **right** of detector center.\n- **Layers** are named for their **radius** from the detector center, like onion layers.\n- **Modules** are named for their **rotation** about the detector cetner.\n\nThis can be seen in the 3D plots above as well.\n\nViewing group distribution:\n"},{"metadata":{"_cell_guid":"17c1b693-4a48-4935-946d-0f2efe4c9fd8","_uuid":"14eabe6b584f2a4353e5925d710e4a368f36f602","trusted":true},"cell_type":"code","source":"hits.volume_id.value_counts()","execution_count":45,"outputs":[]},{"metadata":{"_cell_guid":"de8f4c99-0a56-4361-bd92-2c0de4ea645b","_uuid":"4a2ce8399e4073c322ef03a63f8e591e3b9c3b1d","trusted":true},"cell_type":"code","source":"hits.layer_id.value_counts()","execution_count":46,"outputs":[]},{"metadata":{"_cell_guid":"01c3ea8d-5afb-4168-b194-a6940f8baea5","_uuid":"1a81943ab848f155762d0d5e0f2b0fdf795d35ea","trusted":true},"cell_type":"code","source":"hits.module_id.value_counts().head()","execution_count":47,"outputs":[]},{"metadata":{"collapsed":true,"_cell_guid":"87be901b-9e51-4075-93d6-ded8b6a3093d","_uuid":"cf0098e0387de6e8d8529c4e0d6e7681a0f63968","trusted":false},"cell_type":"markdown","source":"## Hit Feature Correlations\nPlotting each feature (x,y,z, volume, layer, module) against each other."},{"metadata":{"trusted":true,"_uuid":"fb0b4c9f7352417f53823480fdd4905a2eca2f2d"},"cell_type":"code","source":"# Pairplotting 120k hits takes too long - so a random sample of 3k.\nsample = hits.sample(3000)\n# Color coding by group\nsns.pairplot(sample, hue='volume_id', size=8)\nplt.show()","execution_count":48,"outputs":[]},{"metadata":{"_uuid":"fcc52f9f905f3e572b9ec07769d041b29a96a624"},"cell_type":"markdown","source":"Qualitatively, there are correlations between hit features. Usebale ones may be easily extractable.\n\nPlotting correlation heatmap, dropping `hits_id`:"},{"metadata":{"trusted":true,"_uuid":"38c11a49eecea4c7d8e2dcd3e7eb2dab012d603e"},"cell_type":"code","source":"fg,ax = plt.subplots(figsize=(10,10))\nhitscorr = hits.drop('hit_id', axis=1).corr()\nsns.heatmap(hitscorr, cmap='coolwarm', square=True, ax=ax)\nax.set_title('Hits Correlation Heatmap');","execution_count":52,"outputs":[]},{"metadata":{"_uuid":"d7209cc3e1aeb04854ccf6ff6b660b29cce53fcc"},"cell_type":"markdown","source":"`module_id` is correlated with `layer_id`, and `layer_id` is correlated with `volume_id`."},{"metadata":{"_uuid":"0dcf2e26ec8a6d95913daffe16c47705203a7526"},"cell_type":"markdown","source":"## <a name=\"cells\">Cells</a>"},{"metadata":{"trusted":true,"_uuid":"840fc635c585a2a371b9476bf2238cdfb106e1ff"},"cell_type":"code","source":"cells.head()","execution_count":53,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"9433e0630345adfdb2cb4f36fd627fca3ce22f5b"},"cell_type":"code","source":"cells.tail()","execution_count":54,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"f9db6fff7ce81cc550ea148734e8e793483f2631"},"cell_type":"code","source":"cells.describe()","execution_count":55,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"19af63990b558c11c27b7b0ae4a9b02a3b91d32b"},"cell_type":"code","source":"fig,axe = plt.subplots(figsize=(10,10))\ncellscorr = cells.drop('hit_id', axis=1).corr()\nsns.heatmap(cellscorr, cmap='coolwarm', square=True, ax=axe)\naxe.set_title('Cells Correlation Heatmap');","execution_count":57,"outputs":[]},{"metadata":{"_uuid":"fe1e31a2baf71dd33939a9cb95e9cde4066f20bb"},"cell_type":"markdown","source":"## <a name=\"particles\">Particles</a>"},{"metadata":{"trusted":true,"_uuid":"d5651b708f9e228ca3a22abd7238850712311232"},"cell_type":"code","source":"particles.head()","execution_count":58,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"1e1481326f53e7c66f7efb1b52ade314ae16ff14"},"cell_type":"code","source":"particles.tail()","execution_count":59,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"f6f9cfcdfaa3810459e9c480197d07104a9c2552"},"cell_type":"code","source":"particles.describe()","execution_count":60,"outputs":[]},{"metadata":{"_uuid":"091795745f5dc4e5624a239dc5fefdbb515fd75d"},"cell_type":"markdown","source":"Considerations:\n- The particle vertex doesn't necessarily have to be the center of the detector.\n- Charqe `q` is always ±1. No ions.\n- `nhits` can be as low as `0`.\n\nPlotting histograms:"},{"metadata":{"trusted":true,"_uuid":"f8119dfde3363f308fe8b5deab91b020fff24e9a"},"cell_type":"code","source":"plt.figure(figsize=(15,10)); \nplt.subplot(1,2,1); plt.xlabel('Charge (e)'); plt.ylabel('Counts')\nparticles.q.hist(bins=3)\nplt.subplot(1,2,2); plt.xlabel('nhits')\nparticles.nhits.hist(bins=particles.nhits.max())\nplt.show();","execution_count":65,"outputs":[]},{"metadata":{"_uuid":"7b6da1d31ac6efb402716e137151fea4498d93ee"},"cell_type":"markdown","source":"There's about a 20% difference between +q and -q particles. `nhits` seems to be a a distribution around zero + a gaussian around 12. `nhits` and `momentum` should be able to be convolved.\n\nChecking if `nhits` is proportional to total momentum `p`:"},{"metadata":{"trusted":true,"_uuid":"8806d2b4eb17b9717c30af66a6ab2338a45fb557"},"cell_type":"code","source":"fig,axe = plt.subplots(figsize=(10,10))\np = np.sqrt(particles.px**2 + particles.py**2 + particles.pz**2)\naxe.scatter(particles.nhits, p)\naxe.set_yscale('log'); axe.set_xlabel('nhits'); axe.set_ylabel('Momentum [GeV/c]')\nplt.show();","execution_count":68,"outputs":[]},{"metadata":{"_uuid":"36aca3c49e9adb0265b6fe5749f3825380dd5b91"},"cell_type":"markdown","source":"No correlation is apparent. Looking at a histogram of momentum (log axis):"},{"metadata":{"trusted":true,"_uuid":"e3e8b0fa4afbcb6b4e6e3beaae0f8baa954e15c6"},"cell_type":"code","source":"fig,axes = plt.subplots(1,2,figsize=(15,8))\naxes[0].hist(np.sqrt(particles.px**2 + particles.py**2), bins=100, log=True)\naxes[0].set_xlabel('Transverse momentum [GeV/c]'); axes[0].set_ylabel('Counts')\naxes[1].hist(particles.pz.abs(), bins=100, log=True)\naxes[1].set_xlabel('Z momentum [GeV/c]')\nplt.show();","execution_count":69,"outputs":[]},{"metadata":{"_uuid":"0ad798a92c99f4709f07e0b533338b177725fd08"},"cell_type":"markdown","source":"Comparing Transverse and Z momenta:"},{"metadata":{"trusted":true,"_uuid":"adbc5017879e731dac0547c36cf2091f239365cb"},"cell_type":"code","source":"fig,axe = plt.subplots(figsize=(10,10))\naxe.scatter(np.sqrt(particles.px**2 + particles.py**2), particles.pz, s=1)\naxe.set_xscale('log'); axe.set_xlabel('Transverse momentum [GeV/c]'); axe.set_ylabel('Z momentum [GeV/c]')\nplt.show();","execution_count":70,"outputs":[]},{"metadata":{"_uuid":"82ad5fcad5b3d3526cdc308c5c80636d6099929b"},"cell_type":"markdown","source":"Clipping out outliers, such as the pz=500 particle at top right."},{"metadata":{"trusted":true,"_uuid":"4c56833177954e37ab4f247ae5ec8aae07b2b268"},"cell_type":"code","source":"p = particles[particles.pz < 200]\n\nfig,axe = plt.subplots(figsize=(10,10))\naxe.scatter(np.sqrt(p.px**2 + p.py**2), p.pz, s=3, alpha=0.5)\naxe.plot([.1,.1],[p.pz.min(), p.pz.max()], c='g')\naxe.plot([.1,np.sqrt(p.px**2 + p.py**2).max()], [.1,.1], c='r', linestyle='--')\naxe.set_xscale('log'); axe.set_xlabel('Transverse momentum [GeV/c]'); axe.set_ylabel('Z momentum [GeV/c]')\nplt.show();","execution_count":73,"outputs":[]},{"metadata":{"_uuid":"9d10cbd295893bd7f5b2a84adc6fc3598d680486"},"cell_type":"markdown","source":"Plot of momentum in the Z direction (down the length of the accelerator) and transverse (verticle). $\\Rightarrow$ particles on the green line have a perfectly parallel trajectory to the beamline, and on red: perpendicular. This plot shows that  particles spread in a cone-like fashion.\n\nPlotting correlation:"},{"metadata":{"trusted":true,"_uuid":"fc0238e37fee47bb20dc3412632e9d00859e599a"},"cell_type":"code","source":"f,axe = plt.subplots(figsize=(10, 10))\nparticlescorr = particles.drop('particle_id', axis=1).corr()\nsns.heatmap(particlescorr, cmap='coolwarm', square=True, ax = axe)\naxe.set_title('Particles Correlation Heatmap')\nplt.show();","execution_count":75,"outputs":[]},{"metadata":{"_uuid":"91f0eefa205bb451cf83c046f3418f22c953c621"},"cell_type":"markdown","source":"## <a name=\"truth\">Truth</a>\n\nEach entry maps 1 hit to 1 particle."},{"metadata":{"trusted":true,"_uuid":"4cbdd62f65276a2ca7629247f4ec6c291345f53d"},"cell_type":"code","source":"truth.head()","execution_count":76,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"b960be0b525910f9c5ebccc1b72de8c991e8d841"},"cell_type":"code","source":"truth.tail()","execution_count":77,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"dd04958b1fa313d091818132c7a51c17afda1baf"},"cell_type":"code","source":"# looking at a particle\ntruth[truth.particle_id == 22525763437723648]","execution_count":78,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"352ac15f36562a4670f81c8e19d0451e17f253b3"},"cell_type":"code","source":"# Number of unique particles\nlen(truth.particle_id.unique())","execution_count":79,"outputs":[]},{"metadata":{"_uuid":"93004edb66e699048808006bf7b39f570979a3bd"},"cell_type":"markdown","source":"## Plotting Particle Tracks"},{"metadata":{"trusted":true,"_uuid":"95bfb125d58ea4e85c0dd23c5425bdb5c9c10f7a"},"cell_type":"code","source":"# get every kth particle\nk = 100; tracks = truth.particle_id.unique()[1::k]\n\nf,axe = plt.subplots(figsize=(15,15))\nax = f.add_subplot(1,1,1, projection='3d')\nfor track in tracks:\n    t = truth[truth.particle_id == track]\n    ax.plot3D(t.tz, t.tx, t.ty)\nax.set_xlabel('z [mm]'); ax.set_ylabel('x [mm]'); ax.set_zlabel('y [mm]')\n# These two added to widen the 3D space\nax.scatter(3000,3000,3000, s=0); ax.scatter(-3000,-3000,-3000, s=0)\nplt.show();","execution_count":83,"outputs":[]},{"metadata":{"_uuid":"83c548bba774ad55fdec8278e8c3e9b2b2379167"},"cell_type":"markdown","source":"Many particles don't start at the detector cetner, but originate somewhere else. Also more-helical trajectories have less z-momentum."},{"metadata":{"trusted":true,"collapsed":true,"_uuid":"26bd44828c19cfeb666a2d7de2ffa5bfd14dbbb2"},"cell_type":"code","source":"","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"collapsed":true,"_uuid":"fac6f46e3af5b73134f6d38eb4716d188e1add88"},"cell_type":"code","source":"","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"collapsed":true,"_uuid":"15426856c6225865ea464960a21dd998a2d27639"},"cell_type":"code","source":"","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"collapsed":true,"_uuid":"01970dffcce7c57b9aee9769931a6117b9986349"},"cell_type":"code","source":"","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"collapsed":true,"_uuid":"77da942ac5f01c0ff79ac29970fdaf3921bc4765"},"cell_type":"code","source":"","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"collapsed":true,"_uuid":"0fbbb79775f5afe09918bbaa4145d210e2600fb8"},"cell_type":"code","source":"","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"collapsed":true,"_uuid":"6b75284bfc06aa5713b91cd12160d3269ec92b2f"},"cell_type":"code","source":"","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}