{
"cells": [
{
"cell_type": "markdown",
"id": "b1841487",
"metadata": {},
"source": [
"# df016_vecOps\n",
"Process collections in RDataFrame with the help of RVec.\n",
"\n",
"This tutorial shows the potential of the VecOps approach for treating collections\n",
"stored in datasets, a situation very common in HEP data analysis.\n",
"\n",
"\n",
"\n",
"\n",
"**Author:** Danilo Piparo (CERN) \n",
"This notebook tutorial was automatically generated with ROOTBOOK-izer from the macro found in the ROOT repository on Tuesday, March 19, 2024 at 07:06 PM."
]
},
{
"cell_type": "code",
"execution_count": 1,
"id": "c6431ff0",
"metadata": {
"collapsed": false,
"execution": {
"iopub.execute_input": "2024-03-19T19:06:59.546197Z",
"iopub.status.busy": "2024-03-19T19:06:59.545790Z",
"iopub.status.idle": "2024-03-19T19:07:00.668285Z",
"shell.execute_reply": "2024-03-19T19:07:00.667334Z"
}
},
"outputs": [],
"source": [
"using namespace ROOT;"
]
},
{
"cell_type": "markdown",
"id": "84b1ec54",
"metadata": {},
"source": [
"We re-create a set of points in a square.\n",
"This is a technical detail, just to create a dataset to play with!"
]
},
{
"cell_type": "code",
"execution_count": 2,
"id": "fda4a407",
"metadata": {
"collapsed": false,
"execution": {
"iopub.execute_input": "2024-03-19T19:07:00.672949Z",
"iopub.status.busy": "2024-03-19T19:07:00.672618Z",
"iopub.status.idle": "2024-03-19T19:07:03.023908Z",
"shell.execute_reply": "2024-03-19T19:07:03.021647Z"
}
},
"outputs": [],
"source": [
"auto unifGen = [](double) { return gRandom->Uniform(-1.0, 1.0); };\n",
"auto vGen = [&](int len) {\n",
" RVecD v(len);\n",
" std::transform(v.begin(), v.end(), v.begin(), unifGen);\n",
" return v;\n",
"};\n",
"RDataFrame d(1024);\n",
"auto d0 = d.Define(\"len\", []() { return (int)gRandom->Uniform(0, 16); })\n",
" .Define(\"x\", vGen, {\"len\"})\n",
" .Define(\"y\", vGen, {\"len\"});"
]
},
{
"cell_type": "markdown",
"id": "04831334",
"metadata": {},
"source": [
"Now we have in our hands d, a RDataFrame with two columns, x and y, which\n",
"hold collections of coordinates. The sizes of these collections vary.\n",
"Let's now define radii from the x and y coordinates. We'll do it treating\n",
"the collections stored in the columns without looping on the individual elements."
]
},
{
"cell_type": "code",
"execution_count": 3,
"id": "cc3ee114",
"metadata": {
"collapsed": false,
"execution": {
"iopub.execute_input": "2024-03-19T19:07:03.034940Z",
"iopub.status.busy": "2024-03-19T19:07:03.034511Z",
"iopub.status.idle": "2024-03-19T19:07:03.357685Z",
"shell.execute_reply": "2024-03-19T19:07:03.356489Z"
}
},
"outputs": [],
"source": [
"auto d1 = d0.Define(\"r\", \"sqrt(x*x + y*y)\");"
]
},
{
"cell_type": "markdown",
"id": "468030ad",
"metadata": {},
"source": [
"Now we want to plot 2 quarters of a ring with radii .5 and 1.\n",
"Note how the cuts are performed on RVecs, comparing them with integers and\n",
"among themselves."
]
},
{
"cell_type": "code",
"execution_count": 4,
"id": "0d5d6598",
"metadata": {
"collapsed": false,
"execution": {
"iopub.execute_input": "2024-03-19T19:07:03.362177Z",
"iopub.status.busy": "2024-03-19T19:07:03.361835Z",
"iopub.status.idle": "2024-03-19T19:07:07.572891Z",
"shell.execute_reply": "2024-03-19T19:07:07.571554Z"
}
},
"outputs": [],
"source": [
"auto ring_h = d1.Define(\"rInFig\", \"r > .5 && r < 1 && x*y < 0\")\n",
" .Define(\"yFig\", \"y[rInFig]\")\n",
" .Define(\"xFig\", \"x[rInFig]\")\n",
" .Histo2D({\"fig\", \"Two quarters of a ring\", 64, -1.1, 1.1, 64, -1.1, 1.1}, \"xFig\", \"yFig\");\n",
"\n",
"auto cring = new TCanvas();\n",
"ring_h->DrawCopy(\"Colz\");\n",
"\n",
"return 0;"
]
},
{
"cell_type": "markdown",
"id": "9a482c41",
"metadata": {},
"source": [
"Draw all canvases "
]
},
{
"cell_type": "code",
"execution_count": 5,
"id": "35abf533",
"metadata": {
"collapsed": false,
"execution": {
"iopub.execute_input": "2024-03-19T19:07:07.599077Z",
"iopub.status.busy": "2024-03-19T19:07:07.598641Z",
"iopub.status.idle": "2024-03-19T19:07:07.961153Z",
"shell.execute_reply": "2024-03-19T19:07:07.959929Z"
}
},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAArgAAAHYCAIAAAApvgy/AAAABmJLR0QAAAAAAAD5Q7t/AAAgAElEQVR4nO3dX6wc133Y8TNhI5C6oWgKFArBuLQgqIigFGBQ5YmwMjOOoZcyBoNCak0ZhkyyzQMfbKSCX2emj4mK2AXYIKlIq4GltBLQEg774jaaGcvgU1WAfTDkRlAUXgRGQ0WyxPIPXLDbh595fDhzfufuzs7uzsz9fh6E5Zyd2Zm5V/f89sz5nV80m80MAACAzy9t+gQAAMBwESgAAAAVgQIAAFARKAAAABWBAgAAUBEoYPTyPI/mkOf5ps8U86qqKkmSKIqqqlrbJ4r1fBwwIn9v0ycArMke7wPk8pMk2fB5zKGqqjRN1/yheZ7XdR3H8R7/PQHaGFHA6OV5PrufbC/L0t24lzsA6XrHMqYi55ll2Ww2G0VkA0wbIwqYLPqYUVvzj28vx5FAGIECsIDOA/g9jvy7XVr4gPN/6JzvlKkD3pNZ9NJWcUNW93OZ5yNG9HAHWMwMmBztdzuOY3NvTLvxZu9G+/CiLMvG/zhxHM9zJu0dy7LMssw9uPdsZcfGp8j5tw/Y2CvLMvu58lname96XfZU3QNqJ9O4h3PekPD5NJ4fNbTPYZ6fi32bu4tcqb2E9p0MnFL7DXJiu54JMAr8KmOC5O91e7t0BlrP1N7Y+Gccx1mWuV1v+DTcgzd2NIsHCnb37B7tzN3uU95pt0jP1zhg47raH+oeUHZ3D+geP9yvz3MnG2cVOFrjGr3nr+3YuEUzJVCQt4V/7toV7frrAYwFv8qYoMCf6UaT22HYjY14wvYW7nG8Gxu872n0uNrZtgMFb0/c2Oj2xOGjzXlArfdtn3M7CGtzO2arPcwjW+aMObwfEdhx5vwI2uffCBTaZ9u4P4ExKgIFTAZZD9iL7CPnuq7NvS+OdqO8sDkC8p7GZDf5pzRpH+HdsT2UPSf5wtp4BO59GNH+0Da5ujiOGweU02vnR7gH9B48z/OyLAOJFXavxnvkn0VRhE+4zftcYH5z7ts428YNlx9x4z3LnBUwQAQK2Fvki6Pb1cVxLH/oG9GD9KCyUeuPzW5dcnvHzpPd8jxvdEg2FumgEQxZcnqNwzauwl5CkiSNmZW7zvVrz5nodkOSJMnzvLHvQvmf83xu4Odu9N8N5jNiYggUsLe43b+dpi4b5Utt469/YCr7PL1Iv31GVVXSO8pak9qqROETExIKpGnaXsKy/eb2VciX5rqu5QjSbS94NU0dEhTtAo6iw7DEMkhzwB5BeiT2Fvcbs/cP/WD/+idJYr/oyyMD6Z47DyoY36OH+U9mNpvleS6jGqIoirIs13br3JjG3pB+l3Qc4K8BsH4ECthz4jiu67qqKvkCKp1BlmVFUdjV/u334yRJZHv7OPP00H0t42OfMvTVE8tNaI/eL8TeJRnqkAGGWStdcE4LnYm8ebMrLgd+N4Ap4dED9hzp3uxsPtkoHU/7C7r3mX37Ddr29o5z9iveGX/tyYydhxPcRzCN09v1OUL7Pe5kBe0C3ec7jaPNec4u70TRNVvyRwyMBYEC9hz377s7L89ubEy4k2Ci0UPbb7ThT2nvqI2NN9IKvPFKwzLTAmyuQaNXS9O0rutdv9zLgwZvUzhyMq3TlhvSnuQ4j8bJr/9Jgfd3Y/0VrYCVIlDAXtQYSGjwzqWv61pqVUtV63m+0doZf7KjTLtrv036SCnaJNI09SYaFEUh0wLkaLarTtN0nm+xEgHYftp+rmyU6zJzTFywrfaG2EsLz6OUG1IUhc1ZsDdk0aDHvWlyQ+wPxbTSMVbHZsm6vxtr+FxgrTa9kAPQv11/t+331103No5pLbOEc3tBoUb/ahcY9q7M6B5qdv86xN6FlRofEV4yubF2UGMZosCO3rfNc0PaZzvPgksz5aZpCyW5vD9lbQnnec6tcUXaRwAjFc26zjwC9qDOORGVU05J8hfa0xLnOfiKkjKWuS55sei+fV3IcLJU3DORcQX+umIaCBSAddMCBYyOd05oVVXy8IhZjZgG5igAQEeygIR3eubyK1ABA0GgAAAd2emZsjylO6mT4SJMBoECAHQkK1TKDEcZXZB60zx0wJQwRwEAAKgYUQAAACoCBQAAoKIoFACgfyxSOXxzzj0gUAAArARz4IZs/kiORw8AgEmx9dNt4Xgsg0ABADAdeZ5L4TF5TTHP5ZEeCQDoXxRtpn9x19UeTimQAZr/B8SIAgBgImzhb+8S2lKRXEqTm93KxGt7hbfbhx3y+MM9jfYpjWad702VrQQATNhG+hdbjV3Kjss/2032v3Me0NY0t5XH7Xbb1PgUdy9pMr7C8X1deAfzfzqBAgCgf5vqBaXzltduZ+x281KkY9czlLeVZdk4mmxvHE3+aQMIe5BGU3v7psz/A+LRAwBg4uzjAPlnkiT26/6ue9kpDnmeS2fvPVrjAYR9HcdxURT2/e5Tj7E8eiBQAABMXLckyaqqGvGE29m7kiSp69p7kMZUSjnCWEIEQaAAAMACOqdRZFkm8URd1+5shoEjUAAATJybMym0AYDGXo23RVHkHZwoisIde9AeQ8hAgl3jYffzHgYCBQDAxEnfbBdfmrOTdvt1+8ImPbpzFxrHTNPUzkKo69rOnRR1Xc8zQ2I4qPUAAJi+2WwWRZEUOIjjOMuyeSYulGWZpqkti2C7/Mb2LMvc5xFZltmgpN1UFMWIhhMMKzMCAFZhUyszzqn9MKKbqqrcOEAWkJYLbzS137BZrMwIAMAvRFHkfo/XVm9cVGBio7epMZthHFa0ksMGbfqOAgCMGVj/0l6ZcXb/Qoq9nH9gvUX7WctcRY/mP5NBDw11M/DxLgDYCwb4p9gWYkiSZEWVouQjtEoTWtNGzP8DGtwPcnkD/O0EgL2GP8UDN/8PiKwHAMBK2KQAjBqBAgBgJRhRGLL5wziyHgAAgIpAAQAAqAgUAACAijkKAID1aeclLro8Yp7nw0kyXB2bzOlerC0z0VgOsvE22bevLFBGFAAA61PXtdRVsrR3epdAtt3ntOV5LtUiqqqysw5t7co0Te2dkVZ5m7RWVSX7um9byirWe9qsSV4UAIyL9qfYGFOWZXu7bMyyzLbKUobyz/Z/7V6No5Vl6R5kpNy7ZIyRK3JvqbyW6layJcsyWW6ysW/gI+Y8GUYUAACbZ7/+2hf2+7FstMUabGFG2ZLnuf3ObWtA20LPI+XWnIzjWB4iNMpVG2PyPPeOGbj79nAfFgtyhscGU9YELgoAxk77U6x1Q/aF/epclqV8RW4czQS/TNuDjH1QQXiLR0gk4W6RmyaX7L7fvUvtg8x5DuMeUaiqqiiKUYeNALDXNLrwRuucj9Xrui6KQr5qV1VV17UxJo7jKIrkCCuq5rA2dtqBe4tk+CTLMnduR5IkaZqWZbmiSx5roCCTXOwAFABgT5HvytU90pvKC+kdRp0ZIZMZy7J0vwlLSDSbzdxLk8BIrrp9HJk6uuTJjDVQMMYkSRKuEAoAmCrpNe1r6TtlskKe52VZFkWxubNbVlEUjRECGTVpjKBLP9gIieI4li22VOaSJzP66l5RFDXuJiXLAGDjtD/F7RID8jfcfb99HUWRTMfztiZJYp842E5RttjDruDKVs7mN1ryrbgR+sxms/bNtPdN/hm4CXuozLQ3UFj0IGO/CQAwNHxnG7i9Xmaa304M3IGL72hNWzsHvNsfPPa+tsuNy09qTXeOf7LQpxhjPsye0poA7EEjnqMAAABWjUABAACoCBQAAICKQAEAAKhGP5mReYsAAKzO6AMFYFGHz76nNe3b/pnWpOUCBI6mZRyE3dy+7d3+oL5L4INun37au/3opcv68ch6APALPHoAAAAqAgUAAKAiUAAArI+ty2DleT7StZYXUlVVnudauWOttV3aKs/z9g1s79t+W2cECgCA9ZHy0O6WoihsgYapsuFRmqbtaEBrrarKvVdS88Lcq58sG6UUdWNf+zbbugwmMwIA1s12de1uTMpGt780G+frteyuVUd0+1H39QbVdW1z9Nrf/m2rRAzS2q5blOe5rRXp3j17H9I0lSLUtoCkd7BhUQQKmCwtH6FbMsJLb/lL1u7bfk7bpVs9Be20b22rR9u/84DW9NJj/tO+dVU9bXNSbQF6kWWZGyhkWWa/N8tGKSYpBf+klGIcx8apYyRbpNXtg0WaprKvtdbLa6mqSs7fGJPneSMCcFvd6pez2axRRlLKTLvBQWNf+6LfSyZQAACslXSW8pW3KIrZbCaBgtvry3vktS0RLJGBfS1HkNF1t2uU+stSmrmv5/TLCH+hX+jrvg0ypLJ2kiRFUcgLN1CwRy6KoizLhc/4fgQKAIB1k36usdH7KEGCAxk/d6cyhL80yyjF+lfka8cl/X65t48ebHiUZZmMOmRZ1rg/dV3bGGsZBAoAgHWT8QDp59zt7QAiz3P5WrxQcoQdmV/zcwdtACMwbCCjAose344ouNkN9jgyoLL8NEZB1gMAYN3yPJf0B7dndZ/Q2xn+EiVIfz9ncoQcczabuQ/4N8i9LplpKK/bIyhua5udomiMqetadrQPI2zUJS96fOYyzRGF9mRRQWEIABiIOI4bHb/0cPYPuDxct0PrsotM7A8f2T50yLLMJkdsVlmW9rrk3GSiop2E0Wj1kuBJhg3k0uSF7GvDiLquGzmoS/Z90fT6Tjv/BaNzpPiRd7tW/sAYs//KIa1Jy24495ha5uDC6y/oZ+cXKA8R0OGKDp54V9vl1tXHO3yQRisPASyEP8UDN/8PiEcPAABARaAAAABUBAoAAEBFoAAAAFQECgAAQDXN9EhMTGAe/pEdf6JEYK/DZ9VEiTOnXtOazn9wYuFzU5I4wrQ0imsn/SdgjDl82V8ewuhlILolawDYgxhRAAAAKgIFAACgIlAAAKyPrfpoLVTEAetHoAAAWJ/G6sLGmKIo5iziMHa2DOZCre0tbiGoRfftgMmM6O7wWXUOXbe5ch9mT3m3H7j4jrbL7WzhuYT7thc9L2P0lZUDsykfPPa+1nTzo0fVpsUXXe6wuvNdZZIjsB62rqO3b5Nuz90o/7QbZXdvWWpzf9HI9ReQ1EhxByld0S6SGWhN09StFh1FkRSOsmswS4FNWwgqsG83jCgAANbKrYBcVZVbado+mIiiyJaKlte2nqQxJk1TCRTyPG9XAbRFpKT3Xem1zE/67DzPZ7NZewTF25okSePqpLyk3A1bBcqW4a6qqq5rGz9p9REXRaAAAFgr+QYsr91K07afk/7S9vFlWdqu0UYYElLYYKLxEbIlz3M3Ctk495t9+5xtq73MqqoadZvcARK3MGb7yO19OyNQAACsm/1a7G60X4WFbJQXMuHR/SIeGE63g/B1XTceYQxE+9FDQ6C1cWdE+6FDj5ijAABYNxkPSJKk8Y2/3YPK8IMMrc/5rF2eR8gIRH+n3LMeZ06UZZmmaVEUcRyv4pKnGShoD2Yojg4AQyB9eV3X7p/lJEnsIwmZXjCbzSRKkG51oeQIeerf61n3Jnwh4YGQ9jzQJEnsbYyiqPdBlGkGCgQE6/HxK09oTd0SIrQ8hUBqw9FLl7WmDzN1zWPN4bMvaE1Xv/J17/Zj3/22tsv1bTW1QcuhMHoKQ+CWGvOk1qDd7cDPDliPOI4b/aUMMNhvetLNSyKA3UWih10PnmVZURQDyXcQMvfQzquwnb288La25Xmepqk7h9E4wzPhfTubZqAAABgm2827jxjsxvYKAe0t5v5vg1rcMLTnDpK1IUMmdqjDpi96W9uSJInjWGKpOI7tBA534mfvZ06gAACYFDutYdMn0tSOaeaJeNqJD403uI8edv3EDsh6AABMzfKrDMFiRAEAMCnDTIkcL0YUAACAKppegoBd/hp90ZIROlQl6GZr54DWFDiHRx7+iXf78w+plSPe+FRNr9CKJmj1KYwxd//nP9SajvybS1pThzIZgaISmsBNePkLA1rJDuPFn+KBm/8HxIgCAABQESgAAAAVgQIAAFARKAAAABWBAgAAUE1wVipTbQMClRFuXFbrBXTQYeq+0VMYbp9WkxGCFRD8zpx6TWu68Lpa60ErjnDgopo+EBDI49AE0iu0tBSjJ0QEftyUgUAvAn+KbcHoXastSwmDRj1ld5f5S0oOkBTablyg22pfN6pAuRvNvXUj2ktftze6yHoAAAxRFEWy5PBsNkuSxJaA8lb9bYcRskUKQCRJYiskjY6UxzTGpGnaDhRs8Se31EWapnaj3AepECEBh72B9rWtOLUkAgUAwFrZni/PcyndZLs92V5VVaD7T+7J87wsS1uZur2jG2doMcemSDmoPM9ns1m76nRVVVmWVffY7bPZTLbYAQP7Nqk/aTc2hl6WQaAAAFgrWxDZ3B8i2CrJ8lU7iqJ2D9o+lLuj9JTyJdsYk6ap/SD3td0iL7yDGWvgft33dureWpr2Mo0xdV3bg9i7Wte1RAmNIKOzaQYKkWLT5wUAe508F0/TVJ5ByEb3gXpd17PZTL5qz3/Yqqrqupav2rPZTIIAOwfCRhLuLlmWychElm14NdL2XA25FnP/cwRpSpKkKAq5Fhk5kNaiKGxcFUWRPKfo5dHDNItCMZkRAAbLdoF5njem1MkQeudjNvpF90l/OxqQT7cj9pvVOHP3nsgQghs5yZmbewMM8lqeQcgbbAAhgytLhgvTDBSgTcXf2vEXLDDBPIUXn3nTuz1QFCCQCxCY8N8hF+DO8U+0Jq3WQ6CgQ2DCf4fshkCyRiD95PpHjy66y90dNYXh+rb/aFud0lKAJUlwYDt1+brcGGDf9XFDmxyq/b3cPpiwX9DbOnxc77xzFObs3e11uffQvpZ7smSgMM1HDwCAAZLESC3Nz9w/52DX5+sSZ8hQgU25NM5wvTEmy7I0Tb2jFPLQQR5AdLmY5diRDHcsxDujQiYiNPIa5IWbHVoUhRzQHSNxJzF0xogCAGB9yrK0swjNvYkCxpg4juUxhPsGbwdfFIXNdLBj7NLr2660LEt5kee57UFd7qyI9qjGGkjHLxdiz1ZSIZIkcW9C4xbJRtlFnkHYRw822rBH7iUMIlAAAKyPLKLQ3u5O1gvMMwsMM2idvfdo7js3Na2t/bl2i/cmeK/de/L9XhGPHgAAgIpAAQAAqCZYFoFaD0Zf/P/uzgMdjqZlFuy/cmjRXcICB9QcPPGu1rSe6hXf+fy3tF2e/+vf1Zo6XGkgLUWrkRH4oMAPKJCsAcyPP8UDR60HAADQAwIFAACgIlAAAAAqAgUAAKAiUAAAAKoJzkqd2FTbw2ff05o6ZBZ0mGxv9Pn2gcn2AVoJBqPnKXRLbdD20oopmOCE/34TSQLOnHrNuz1QouLWVbWER4fTCxS8AOYX+FNsl1tuV2dokCWK3bUFbakI9w09nO7eM39fycqMAID1keoM0tm71SO9/Va7oJHdUV6naWpXcR6dqqrkAr2xjhsPNSpBuBvN/ctRt4+wfCDFowcAwFrZLi3Pc6nm0KgCJUUmtd2Te/I8L8vS1n1o7+j2rI1eNvzPNZAoxxiTpmm7L3cLZNsrStPUbrQnLIWk3ZJRjeMvb5qBQqTY9HkBAH5eFVpeuyGC/QYsPVwURbvWgHa/akvoYO71neb+Gozua3P/WIV7Pmsj9Z/yPJ/NZt4y0zLuIuz22WwmW+wogjy+kSrbbpDUV5RgphoozBSbPi8A2OvkT3GaplEUeQfP67qezWbSg85/2KqqpIC17CjdpJ0D0R6Ez/Pcds91XW/k4YV7Pt5IpV3mSuZk2I1ucW23orTU0vTW3uxgmoECAGCwqqqSctLGmPaAebfuzYYC7iN/qTEtL7Is8+7V+RN71J7UKUGPuVcz2txfXbMoCrnGJEmk9nQURXEcu+MrPYY+TGbcWwKJEls7B7SmD7OnvNu1RAATnmx/Qs166JCn8OOvfF1r+kf1H3i3b6lnZo5euqw3+jMLAiUYHjz2vtb0/EPvaE0XXn9BPwe/QJ5C4GcErJ9847eduvSCje/Nuz5uaJNDtbtb23HaftdlZ0EOIW+icQ7uaIoMIbhDLDIJ1NyLBmS7fVuappSZBgCMkgyVa7P3zf1zDnadNyBxhgwVuIPw7sy+LMvSNPWOGcjTh009d3B55yjMs2NRFO7M0KIoZEcZY6jr2jtTclEECgCA9SnLUiYoiCzLpCeTIXT3DVoHXxSF7OvmRspTebtdnmsYPXXQ2tRzBzv30J0/4Z16KZMP3OjHNrkTGGViY5IkdlpeHMdlWS4/T5NHDwCA9ZGerL3dfQYfGDkPdHvtqX8iPA6/qeEE6fhlCoUNayTESZJEoiXZ2I6l7C5yELtxRXP2CRQAAHvRxicotPt1u8UbLXmDpG5B1UJ49AAA2KNIm5/HpMoiiInVegjMww+s8K/Nt+9WFECb2B+Y1R/IU3jjc3+sNX3th9/o69yMfnqB8hCBoy36KcaYaydPaE0HLqpZD1otjMAtDdCOFrgJ1HpALyb2p3h65v8BMaIAAABUzFEAAKwEC+dPw7ACBbuotVZNy32PCLwTALApe/a5gzukf/W/nlrdB/36s3+2nps8oEBBSl1JSmtRFFrl0Eahi0bcAAAAejSgQEGiBLtItbvgVMOeDVQBAFizAc1KjaJIFprw/lPIqEP4nEc61VZbk//m9u0OR9MKNwTSBwJz3bUp+vuvHFr0xEyw3oR2QK0GRJg2sf/Mqde0XV7+gqdsjDh89r1Fj/bq289pTVr5jMAHBX5A2i4BpDYAKzK9Rw9DyXpo1wA1vsUi7NtkLaqNL9ANAMC0DSVQ8NImH8jyllmW2VKbDdGCVnoVAACM14DmKLS1gwB3KW9plYWyG8b46AEAgAEa9IjCrkiMBABgpYYyomArbLYLkzfeliTJeKcmBCed+acZbuk7BBYPvnFl4Xl8Ry89rTVt7fgXfr6pT0u8fVo9mjZt0xhzV9keWL04MKFSW4/5/AfqysoX9B+QdrTAjMWADtMPA5iZCGB1BjSiEMexXSNBQgQbPbiPG4qicCt2b6qUOAAAe8GAAoWqqqTYdhRFdV278wxscJDnucQT8ja77gIAAFiFoTx6EN5ev1GWW97TeEgBAABWYUAjCgshSgAAYA3GGigAAIA1GNajh2kITGjvdxXeayfV2fvmpH/zgYv6Lrr9ysLPgYyMQGrDi8+8qTW98amSK6EsxhymrX597rHL2i7x5/5Ga/raD7/h3R5YjBkAJoBAAQCAlbj2S7+y6VPoAY8eAACAikABAACoCBQAAICKQAEAAKiYzLhWBy6+ozVtbR/wbg/kCBw++4LWdEepwhCY8H/hdfVoWrJGIFNDq4xggqUWtnb8NyFwtEDBiztKhYjAlb5x4l2t6a6S+gEA08aIAgAAUE1zRCGKIu92dyloAACwq2mOKMwUmz4vAABWKM/zJElsyeVeTDNQAABgr4miqCgKY0xRFD1WRCJQAABg9JIkieN4NptVVVWWZV3X3oLMHUxzjsJmBQo6BCogaCUDjhTPdfigl94qtKYOtGSN/eaQtksgGeF/fOaPtKZfNf/K/0FX1A8yx9QWLYfC6DkUN/SiEoG7DQAbV9d1WZbyOkmSHp+2EygAADBuMnjgzk7ocZoCgQIAABv2pS/+u+UPEkVRHMfGmLqui6Loa1CBQAEAgA373n/754vu0o4tsiyzAwlRFOV53su4ApMZAQCYAjcsiOO4r8mMBAoAAIybJEO6kUFd131lSPLooaNAmYOAfdsLH/DMqS61HvZt+3Mlbm7f1nbZ0nMBtAoR541atcHo6QPHzLe1pkeUUgvXj6ufs+/q42qbQkswAYDxiuM4TVOZlyBDCwQKAADg56qqSpLEVjAoy5JAAQCAQbv2S7+yzo/ra1JCA3MUAACAikABAACoCBQAAICKOQodBVb+1yojmED1Ad35D9TMAi1HIODmR49qTYHqDK++7c+heETfJZCnELgJWqmF/erBzM3jn3T4IADAnBhRAAAAqmmOKNj8kIYeq2kBALAXTDNQICAAAKAXPHoAAAAqAgUAAKCa5qOHHmklGPbplRHe+Pwfa01fvfJNremgksIQqGVw7bSaEKGddodEid4FSi0cKX7k3R4oUbH/yiH1g/TMFADAnBhRAAAAKgIFAACgIlAAAAAqAgUAAKAiUAAAACqyHnahZTfc3XlA2+XcjS8tejRjzHWtCoM+4V9LbTB6DsUtPYciUOuhg9unn9aatNSGbgJFNwAAy2NEAQAAqAgUAACAikcPAACsxF/N1BXhRoQRBQAAoCJQAAAAKh49dBTIX+iWWfDhSX/hhkBqQ2DC/5HCf3qBOgsHLqrpFVvKdjVTY5fTXrjWw7nHLmu7GKOmVwAAljfNQCGKIu/22Wy25jMBAGDUphkoEBAAANAL5igAAAAVgQIAAFARKAAAANU05yj0SKvpEMh6CAikCWjuHP9Ea+pQNOHAxXe0pn6rMwTuTyAhwhj/3X75C9miJwAA6AUjCgAAQEWgAAAAVDx6AABgCvI8D/yzMwIFAABGL8/zoigaW3o5MoHCLrRJeTe31QWPH3n4J1rTTX0y40tvFd7t5x5Tz+3l0+oUP23S4tbOAW2XwIxFbUbn//md59WjXbmkNQXmOWofBAAIq6oqjuOqqno/MnMUAACYgiRJVnFYAgUAAEavruuqqqIoiqIoSZIehxZ49AAAwIb969/KezlOWZbGmDzP0zTtq+wRgQIAABv2L/8iX3SXRmzhhgUytJDneS/zGQkUAABYib+6e3hTH93jxEYChV0Eshs0106e0JqOXrqsNb3xqX8F5ecfUhddDth/5ZB3+019QeiA29lT3u0HLv57dZdXntCaQstIKx8EAAioqirPczcyqOs6y/pZ/H6akxkjxabPCwCA/iVJUte1zXqQJw6soxDS1wwOAABGoSzLNE3tV2KZ1diLaQYKAADsKUmSrOhL8jQfPQAAgF4QKAAAABWPHnahFW64rldtCBRN+M7nv6c11dFnvdsvvP6Cto7If/kAABOnSURBVMv5D9T0ga1tf02HQK2HAO2K9u/4cyuMMUd21JsQ2MucXuS0AACrx4gCAABQESgAAAAVgQIAAFARKAAAABWBAgAAUJH1sAutcEMgteHFZ97Ums59+qVFT+COXp0hkMLwoVqdoUvliNun/XUoAjehG+2A2uUAAFZtlIFCVVVS+iJJEru0NQAA6N34Hj1UVZWmqcQKaZr2VfQCAAC0jS9QSNNUymxXVZVlWVEUmz4jAAAma3yBgnFKZ8oLtwI3AADo0cgCBTs1ob0RAAD0bpSTGRvagYItyD2nDqU5g/PwM63h5UU/pm9a/kI3JCMAS9ISkbQqMyZYaGb/FbWQyr7tn3m3B7K0zn/gz/kKf1AHZ0695t0eqHRz8MS7WtOtq48vegJ3dx7Qmj5+5YlFj+b6y7sPL7P7QEwhUGgnPqyoJjcAAHvNyB49AACAdRpZoCCDB41nDSylAADAiowsUDDGxHGcpqm8lhCBQAEAgBUZ3xyFqqqSJLHTFZmOAADA6owvUDDkQwKYCi194Na2Wsnljc9/S2v62s43tCZtYn8gs+C2PuH/yM7CZVkOn31Pa3r17ee82wOVbg5qDcbc3L6tNWlpX4FzgxnjowcAALA2BAoAAEBFoAAAAFQECgAAQEWgAAAAVKPMegCAadBKMAT89m/9idb04I3LWtONnScX/aAjhT+1wQSLI2gCKQxbO/4UD227McZ0KqGgZTcsWdBh8hhRAAAAKgIFAACgIlAAAAAqAgUAACal3xJI0wwUIsWmzwsAgNVKkqSu6x5rHUwz64FKUQCmKlCY4M7xR7Wmc6de824P1HoIpDZoaQJHL6lpF7dPn9CatPSKB4+9r+1y6+rjWpPRaz3sheyGqqrquu73mNMMFAAA2Lif3vnMmj8xTdMsy4qi6PGY03z0AADAXpMkSZZleZ73e1hGFAAA2LCd3/ntJY+Q53m/UxMsAgUAADZs+z//+aK7uLFFVVVFUaxofh6BAgAA4yaPG9ysyDRN4zjuZYCBQAEAViuQp7Bv27/9w+ypDh/00lvqFLZAdoPmjJIoYYwxJvNuvf6RmnbR4SYEjrZfz8jY0hr09Ipud3tQ8jx3Y4K6ruM47muyAoECAADjliSJO5xQFEWe530tu0TWAwAAUDGiAADApPQ7q5ERBQAAoGJEAQBWK7BysDb98MBFdcHj/VcOaU37tp9b6MTCzn+gnsMFZWbifqOeW2Bq5BufPu3d/p2D31NP7nNqy/N//bta0yMP/0RpGf1kxpViRAEAAKgIFAAAgIpAAQAAqAgUAACAapqTGaMo8m5f0TrYAABM1TQDBQICAMNx9NJlrenaSf9ayEZfjPnVHTW14a6+sPGd4594twdyKAL2bf/Mu/3BY+9ru3RYQ/qcmnVhbl19XGu6ra/HrC4jfXKBs9qDePQAAABUBAoAAEA1zUcPAABs3E/vdHmyMzSMKAAAABWBAgAAUEXTSxCIogleFICBU2fU6zkCAYH8hYDAB2kHPHjiXW2XQGaBdrRAQYcALSGiw+WE9/pQT4jokdsHPfQf3l7dB9348m+up7NjRAEAAKgIFAAAgIpAAQAAqAgUAACAikABAACoWHAJAFarQwrDn37597Wmcze+pDVd/+hRrekRpQpDILXh5vZtrWlL2X7+A7U8Q6CoxMevPOHdfuDiO9ougYIO6B0jCgAAQEWgAAAAVAQKAABANc05ClEUebezYiMAAAuZZqBAQAAAQC+mGSgAwIpoU/H3G3VWf4eaBXX0WW2XG5ef1Jq29A+6tePPbghkZGipDYG99ncqUaEJJEocffiy1nTtpJp5gW6YowAAAFQECgAAQEWgAAAAVMxRAABgJf7vrQc3fQo9YEQBAIApqKoqz/MkSaqq6vGwjCgAQNPhs+9pTbdfedq7/cjOj7RdApkFV7/yde/2Y9/9trbLwRPvak2BhIg7xz/xbu+WpxDI49AEKkdodzvwKYESFebkIqc1IXmeF0URx7ExJk3TOI77ChcYUQAAYPSKoijLsqqqqqqyLKvruq8jEygAADBuMniQJIn8M89zu3F5BAoAAIxbkiR2SeKqqiRisHHDkpijAADAht058xu9HCdJEnnoUJZlLwc0BAoAAGzc/gv/fdFdvLGFPG7I8zxN07IsexlUIFAAsEcdKdQ8hTOn3tSaDp99wbv9znF1Vv/t7Cn9HP5AOZo/ScEYc+ejR7Wm/VqDMece8xdHeHXnOX2nhQXyO4ye9RC4WM3WzoFFd5k8+8TB3MuAcLcsgzkKAACMmwwhrOjgBAoAAIybpDnIf+0LJjMCAICfK8syTdOiKOw/CRRCoijybrfZIwAATImbIdmvaQYKBAQAAPRimoECAOwqMEX/wuv+1AbTqcxBgFYBoeOsfv3cXn3bn93woZ6R8ed/8S+0pq/+2TcXOq+w/VcOebcHbnXgtNE7JjMCAAAVgQIAAFARKAAAABWBAgAAUDGZEQCAlXjgkxV2sndWd+j7ESgAGJbDZ9/TmgLT4LX0AWPMIw//xLv94An1HK6dVNuOXvIXTXhQPZg5cPEdrUk7t+tGL+ig5AiEafUUAnf741f+RGva90O1TIYmkMfx4StPLHo0rBOPHgAAgIpAAQAAqAgUAACAalhzFKqqqqrKGJMkiVbNwr5HBN4JAACWNKBAoaqqNE3jODbGFEWRZZmtmOlqlNxuxA0AAKBH0XDqJ0VRFMex9Pp5nhdF4T23KNrlnHd9A4B+HSnUOfCBZASNlghgjLl19fFFj2b0ugDd0is0gcoR/ZaHWNs5BOopaKkfgR/Qi8+8qTV1qEMxZG4fdOgPr67ugz79vV9fT2c3rDkKdghBXrSHChg8AABgnYYSKNipCe2N3rdFURRFkffZBAAA6MtQAgUvbfwgSZKyLLMsK4rCO5MxWtBKrwIAgPFa62RGbeKhNjDQDgLyPLdvltaiKNo7MkcBAIBeDCjroYMkSbyBAoDOOkzxC8yhM/pkxtunn/ZuD0yNfPDY+1rTjctPak3aFQXm93X4oG7LS2vOPeafLWiMOf+BvvT04h+k/RRM8DfhznH/CtNb+gddeP2FBU4LQ7LWQCGw5oFsr6rKfUP7zXIEpiYAALAeA5qjEMexXSNBQgQbPbiPG4qisM8v7LoLAABgFQYUKFRVFcexzC6s69qdZ2CDgzzPJZ6Qt9l1FwAAwCoMa46Ct9dPksQNGuQ9jYcUAABgFQY0orAQogQAANZgWCMKADauw1q/Z069pjUF5rofMO94t2+ZA9ouzz/k38UYc16Zh2+M2X/lkNakCSxFvJ7Uj8Nn1XPe0n9AgTWPtRSGow+r6RVnTql3W1t0+Tuf/5a2Sx19Vmt6+QuZ1jR2v/zpL2/6FHow1hEFAACwBgQKAABARaAAAABUBAoAAEBFoAAAAFTR9OonRdEELwpYmwMX1bnumq0dNU8hkAtw8MS73u2BjINA0YRAaoOWp/DiM29quwSSNTpkhQRo5xBIBLhW/WOt6dh3v601aZkp3a5U+7F+/MoT2i57h9sHBQqXLO/v8l9rdHbVPf3WOmBEAQCA0cvzPE1TWZOwKIooivo6MoECAACjVxRFlmUyoiAjDX0NKkxzwSUtkuKRBABgqtzIoMdaSNMcUZgpNn1eAACsRKOPq+u6r1oH0xxRAABgRP4u/7W+DlVVVZqmpr9HDxNMECDrAVhGIOtBy27oloygzZA/ekmtPhCo9RBIE3jprcK7XStYYILJGv1mPWjVGQI/Ba08hAnOsQ+VolAErlT7iXf4cU/PBrMejDFJktR13eNzB8OIAgAA0yDz88qy7LfAMoECAACjF0VRvwMJFoECAADjJvFBkiSNQKGXoQUCBQAAxs2us1QUv5iR09cAwzTTIwEA2DvyPG+vCNDXYwhGFADcJzBx/a62y+Iz6k1oQrha6+H89qPq4ZTUBmPM+Q9OeLd3O20tfSCQI6ClNhg9xWP/lSe1XY7sdJlIryUdBNIrAr8JRsl62DupDXsKgQIAACvx/z7Zt+lT6AGPHgAAgIpAAQAAqAgUAACAikABAAComMwI4D6B2fsvPvOmd/uF11/o8EFa+kBg5nygDEToHI5/4t38p1/+fW2Pcze+pDVdO+nPoQis6n/47Hta05lT/qSDC0bNeghU1tCKcZjA6W1re5DCgJ9jRAEAAKimOaIghTHaqCoJAMBCphkoEBAAANALHj0AAAAVgQIAAFBF0xulj6IJXhSwNoEp+gdPvOvdfv0jtQRDoF6All6hZUOY4Dz8wGmvR7dkjRuX/dkN2q02xty6qtbCePDY+1pTh2SNQIkKBLh90MO/9+PVfdDHf/jkejo7RhQAAICKQAEAAKgIFAAAgIpAAQAAqAgUAACAapoLLgFYhR9+5o+8249d/ra2S6ByhJbdcObUa9ouRy89rX+Qmgugzd4/cNFfZyHs9mn/OQTSB4xRz027P1o2hAkmRARoWSEfv0JqA3bBiAIAAFAxogAAwErMPt236VPoASMKAABARaAAAABUBAoAAEDFHAUA8zqa/Bfv9jvvq+kDWzsHFv2Ul7+QaU2v6pkFgQoRWqmF26f95Q9MMCEimN3gFyiaoH3Qln60QEJEoN6EOakfEQiaZqAQRZF3O8WiAABYyDQDBQICAAB6wRwFAACgIlAAAAAqAgUAAKCa5hwFAJ0FqjNoE/5vZ2oJhkD6wH4lT0GrSmCMuXP8ttZ07pk3taZX337Ou/2lhwptF2PUhAgtveLO8U+0XUI34coh/6foZxZKbQCMyfM8z/MeD8iIAgAAE1FVVVEUVVX1eEwCBQAARq+qqiRJ0jTt/cgECgAATEGSJFmmrlfWGXMUAAAYvSRJkiQxxhRFYPJNFwQKAO4TWG/4pbe0P0Dqlxhtsp4x5n//W/8UyEPffUvb5fZpddakMWrT+Q/8cwnf+FTd5ZGHf6I13dr2L0oduNLQ9MPTagv2lJ9e+AebPgUVgQIAABv2mTN/uegua4stmKMAAABUBAoAAEDFowcAAFbigY/VZbhGhBEFAACgYkQBwH0C6w2/8bm/8W7XlnY2xhw88b7WpGU3BNIHDl9RV3cOrDxttv2bn39IvdJ45r9SY8xXd7658AkAazSbzfo94DQDhSiKvNt7v30AAEzbNAMFAgIAAHrBHAUAAKAiUAAAACoCBQAAoIqm9zg/iiZ4UcDaHD67cGZBoDzE0UuXtaZbVx/3br+784C2y8ET72pNNy4/qTVppx34oFB1BiDI7YP+/j9RM2uW97f/6TfW09kxogAAAFQECgAAQEWgAAAAVAQKAABARaAAAABU01yZEUBnZ069pjWd/+CEd3ug1oMx/tQGo+dKBIpNfHzSfwLGGHNSPwUAS2BEAQAAqEYfKOR5vulTAABgssYdKFRVVRRFVVWbPhEAAKZprIFCVVVJkqRpuukTAQBgykY8mTFJkiRJiqLY9IkAAODxmevqgujL+9vVHfp+oy+LEEVRWZZJkrhbxn5RAICRcvugX/3N/7i6D/pfb/+z9XR2Ix5RCIiiaKH3E1gAAOA13EChqirvLMV50hzo+AEA6MVYJzMCAIA1GO6IgsxV3PRZAACwpzGiAAAAVAQKIYtOikRfuPMbxM3fIG7+pnDnA4b76GFOzFsEAGB1GFEAAAAqAgUAAKAiUAAAACoCBQAAoCJQWJN+p9T2eLRhHqpfw7zGwf5K9Giw1zjMQ/VrmNc42F8JBIw+6wEAABin9EG/KxYyogAAwOhVVZWmqcQKaZrOUxdpTgQKAACMXpqmcRxLoJBlWVEUfR2ZQAEAgCmwowjywluBuQMCBQAAxs1OTWhvXN40JzMOc47uYKf7Tv5Q/R5tmIfq92jDPFS/Rxvmofo92uQP1e/RVpFD8eMf/NMNfjqBgorqDwCAcVlFz9VX4gOPHgAAgIpAAQCAcZPBg8azhr5GFCIG6ueX53mPmalwzbNOiH2P6HdFEZiVrdaCBn7bh2nsf+GTJKnrWvp093UPZphPWZbGmLIsN30iEyT3No7jOI6NMVmWed/W+NWN43itZzl1c/4UsCR+24dpGn/h5Zeq986dQGF3ZVnauz/2X6Nhcv8OZlmm/YoT167UnD8FLInf9qHhL/yumKMwlyRJ5H9prMiu64T0leeDgBWt1oIGftuHhr/wYQQKu0uSZOzProZsznVC7NuiKIqiiB9Hv1a6WgssftsHiL/wuyJQwBBpXVSSJGVZyjLmzO1aNQKF9eC3HQM3wQWXOmhMMLaIMXv0/e9//wc/+EF7+xe/+MX2xvafRTfkl9YeS57Ai85pPfhtx8AxooBRog/D3sFvOzaLEQVjyFFei2efffbZZ5/VWquqcn8E7R+H/IwY41kRu1pL+KeAJc15n/ltx6AwooDNi+M4TVN5LX807d9TdwC2KAr7hEgqr6/7RCdN+ymgX/y2Y3w2nZ85JoYs25XxrhPSWALFfQ/rz6zCilZrQQO/7cPEX3gNSzhjZBrDtsCE8duOISBQAAAAKuYoAAAAFYECAABQESgAAAAVgQIAAFARKAAAABWBAgAAUBEoAAAAFYECAABQESgAAAAVgQIAAFARKAAAABWBAgAAUBEoAAAAFYECAABQESgAAAAVgQIAAFARKAAAABWBAgAAUBEoAAAAFYECAABQESgAAAAVgQIAAFARKAAAABWBAgAAUBEoAAAAFYECAABQESgAAAAVgQIAAFARKAAAABWBAgAAUP1/6aFy0lcT9y0AAAAASUVORK5CYII=\n",
"text/plain": [
""
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"gROOT->GetListOfCanvases()->Draw()"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "ROOT C++",
"language": "c++",
"name": "root"
},
"language_info": {
"codemirror_mode": "text/x-c++src",
"file_extension": ".C",
"mimetype": " text/x-c++src",
"name": "c++"
}
},
"nbformat": 4,
"nbformat_minor": 5
}