{ "cells": [ { "cell_type": "markdown", "id": "ef458ecb", "metadata": {}, "source": [ "# hist102_TH2_contour_list\n", "\n", "#### Image produced by `.x ContourList.C`\n", "The contours values are drawn next to each contour.\n", "\n", "#### Output produced by `.x ContourList.C`\n", "It shows that 6 contours and 12 graphs were found.\n", "\n", "#### `ContourList.C`\n", "\n", "Olivier Couet\n", "\n", "\n", "**Author:** Josh de Bever (CSI Medical Physics Group, The University of Western Ontario, London, Ontario, Canada), \n", "This notebook tutorial was automatically generated with ROOTBOOK-izer from the macro found in the ROOT repository on Tuesday, May 19, 2026 at 08:13 PM." ] }, { "cell_type": "code", "execution_count": null, "id": "9dfb4946", "metadata": { "collapsed": false }, "outputs": [], "source": [ "Double_t SawTooth(Double_t x, Double_t WaveLen);" ] }, { "cell_type": "markdown", "id": "f2f4e865", "metadata": {}, "source": [ " Definition of a helper function: " ] }, { "cell_type": "code", "execution_count": null, "id": "738032e5", "metadata": { "collapsed": false }, "outputs": [], "source": [ "%%cpp -d\n", "Double_t SawTooth(Double_t x, Double_t WaveLen)\n", "{\n", "\n", " // This function is specific to a sawtooth function with period\n", " // WaveLen, symmetric about x = 0, and with amplitude = 1. Each segment\n", " // is 1/4 of the wavelength.\n", " //\n", " // |\n", " // /\\ |\n", " // / \\ |\n", " // / \\ |\n", " // / \\|\n", " // /--------\\--------/------------\n", " // |\\ /\n", " // | \\ /\n", " // | \\ /\n", " // | \\/\n", " //\n", "\n", " Double_t y;\n", " if ((x < -WaveLen / 2) || (x > WaveLen / 2))\n", " y = -99999999; // Error X out of bounds\n", " if (x <= -WaveLen / 4) {\n", " y = x + 2.0;\n", " } else if ((x > -WaveLen / 4) && (x <= WaveLen / 4)) {\n", " y = -x;\n", " } else if ((x > WaveLen / 4) && (x <= WaveLen / 2)) {\n", " y = x - 2.0;\n", " }\n", " return y;\n", "}" ] }, { "cell_type": "code", "execution_count": null, "id": "6211e285", "metadata": { "collapsed": false }, "outputs": [], "source": [ "const Double_t PI = TMath::Pi();\n", "\n", "TCanvas *c = new TCanvas(\"c\", \"Contour List\", 0, 0, 600, 600);\n", "c->SetRightMargin(0.15);\n", "c->SetTopMargin(0.15);\n", "\n", "Int_t i, j;\n", "\n", "Int_t nZsamples = 80;\n", "Int_t nPhiSamples = 80;\n", "\n", "Double_t HofZwavelength = 4.0; // 4 meters\n", "Double_t dZ = HofZwavelength / (Double_t)(nZsamples - 1);\n", "Double_t dPhi = 2 * PI / (Double_t)(nPhiSamples - 1);\n", "\n", "TArrayD z(nZsamples);\n", "TArrayD HofZ(nZsamples);\n", "TArrayD phi(nPhiSamples);\n", "TArrayD FofPhi(nPhiSamples);" ] }, { "cell_type": "markdown", "id": "8e274b48", "metadata": {}, "source": [ "Discretized Z and Phi Values" ] }, { "cell_type": "code", "execution_count": null, "id": "8c968984", "metadata": { "collapsed": false }, "outputs": [], "source": [ "for (i = 0; i < nZsamples; i++) {\n", " z[i] = (i)*dZ - HofZwavelength / 2.0;\n", " HofZ[i] = SawTooth(z[i], HofZwavelength);\n", "}\n", "\n", "for (Int_t i = 0; i < nPhiSamples; i++) {\n", " phi[i] = (i)*dPhi;\n", " FofPhi[i] = sin(phi[i]);\n", "}" ] }, { "cell_type": "markdown", "id": "5b1c6dcf", "metadata": {}, "source": [ "Create Histogram" ] }, { "cell_type": "code", "execution_count": null, "id": "66f39906", "metadata": { "collapsed": false }, "outputs": [], "source": [ "TH2D *HistStreamFn =\n", " new TH2D(\"HstreamFn\",\n", " \"#splitline{Histogram with negative and positive contents. Six contours are defined.}{It is plotted \"\n", " \"with options CONT LIST to retrieve the contours points in TGraphs}\",\n", " nZsamples, z[0], z[nZsamples - 1], nPhiSamples, phi[0], phi[nPhiSamples - 1]);" ] }, { "cell_type": "markdown", "id": "8d7c2ffc", "metadata": {}, "source": [ "Load Histogram Data" ] }, { "cell_type": "code", "execution_count": null, "id": "c0c2d15b", "metadata": { "collapsed": false }, "outputs": [], "source": [ "for (Int_t i = 0; i < nZsamples; i++) {\n", " for (Int_t j = 0; j < nPhiSamples; j++) {\n", " HistStreamFn->SetBinContent(i, j, HofZ[i] * FofPhi[j]);\n", " }\n", "}\n", "\n", "gStyle->SetOptStat(0);\n", "gStyle->SetTitleW(0.99);\n", "gStyle->SetTitleH(0.08);\n", "\n", "Double_t contours[6];\n", "contours[0] = -0.7;\n", "contours[1] = -0.5;\n", "contours[2] = -0.1;\n", "contours[3] = 0.1;\n", "contours[4] = 0.4;\n", "contours[5] = 0.8;\n", "\n", "HistStreamFn->SetContour(6, contours);" ] }, { "cell_type": "markdown", "id": "6b171287", "metadata": {}, "source": [ "Draw contours as filled regions, and Save points" ] }, { "cell_type": "code", "execution_count": null, "id": "b5a5500c", "metadata": { "collapsed": false }, "outputs": [], "source": [ "HistStreamFn->Draw(\"CONT Z LIST\");\n", "c->Update(); // Needed to force the plotting and retrieve the contours in TGraphs" ] }, { "cell_type": "markdown", "id": "1625d62c", "metadata": {}, "source": [ "Get Contours" ] }, { "cell_type": "code", "execution_count": null, "id": "8643f34c", "metadata": { "collapsed": false }, "outputs": [], "source": [ "TObjArray *conts = (TObjArray *)gROOT->GetListOfSpecials()->FindObject(\"contours\");\n", "\n", "if (!conts) {\n", " printf(\"*** No Contours Were Extracted!\\n\");\n", " return;\n", "}\n", "\n", "TList *contLevel = nullptr;\n", "TGraph *curv = nullptr;\n", "TGraph *gc = nullptr;\n", "\n", "Int_t nGraphs = 0;\n", "Int_t TotalConts = conts->GetSize();\n", "\n", "printf(\"TotalConts = %d\\n\", TotalConts);\n", "\n", "for (i = 0; i < TotalConts; i++) {\n", " contLevel = (TList *)conts->At(i);\n", " printf(\"Contour %d has %d Graphs\\n\", i, contLevel->GetSize());\n", " nGraphs += contLevel->GetSize();\n", "}\n", "\n", "nGraphs = 0;\n", "\n", "TCanvas *c1 = new TCanvas(\"c1\", \"Contour List\", 610, 0, 600, 600);\n", "c1->SetTopMargin(0.15);\n", "TH2F *hr = new TH2F(\n", " \"hr\",\n", " \"#splitline{Negative contours are returned first (highest to lowest). Positive contours are returned \"\n", " \"from}{lowest to highest. On this plot Negative contours are drawn in red and positive contours in blue.}\",\n", " 2, -2, 2, 2, 0, 6.5);\n", "\n", "hr->Draw();\n", "Double_t xval0, yval0, zval0;\n", "TLatex l;\n", "l.SetTextSize(0.03);\n", "char val[20];\n", "\n", "for (i = 0; i < TotalConts; i++) {\n", " contLevel = (TList *)conts->At(i);\n", " if (i < 3)\n", " zval0 = contours[2 - i];\n", " else\n", " zval0 = contours[i];\n", " printf(\"Z-Level Passed in as: Z = %f\\n\", zval0);\n", "\n", " // Get first graph from list on curves on this level\n", " curv = (TGraph *)contLevel->First();\n", " for (j = 0; j < contLevel->GetSize(); j++) {\n", " curv->GetPoint(0, xval0, yval0);\n", " if (zval0 < 0)\n", " curv->SetLineColor(kRed);\n", " if (zval0 > 0)\n", " curv->SetLineColor(kBlue);\n", " nGraphs++;\n", " printf(\"\\tGraph: %d -- %d Elements\\n\", nGraphs, curv->GetN());\n", "\n", " // Draw clones of the graphs to avoid deletions in case the 1st\n", " // pad is redrawn.\n", " gc = (TGraph *)curv->Clone();\n", " gc->Draw(\"C\");\n", "\n", " l.DrawLatex(xval0, yval0, Form(\"%g\", zval0));\n", " curv = (TGraph *)contLevel->After(curv); // Get Next graph\n", " }\n", "}\n", "c1->Update();\n", "printf(\"\\n\\n\\tExtracted %d Contours and %d Graphs \\n\", TotalConts, nGraphs);\n", "gStyle->SetTitleW(0.);\n", "gStyle->SetTitleH(0.);" ] }, { "cell_type": "markdown", "id": "b353ca4a", "metadata": {}, "source": [ "Draw all canvases " ] }, { "cell_type": "code", "execution_count": null, "id": "3a7f0126", "metadata": { "collapsed": false }, "outputs": [], "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 }