01-introduction.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# PartiQL query language demonstrator\n",
    "\n",
    "In our discussions about domain-specific languages for particle physics, we agreed that we want the language to be \"declarative,\" to \"specify *what* is to be computed, rather than *how* it is to be computed.\" However, we haven't talked much about specific language features that would achieve this end.\n",
    "\n",
    "PartiQL is a toy language intended to demonstrate several features that would be a radical departure from general-purpose languages, addressing problems specific to particle physics. My intent is to inject ideas into the development of languages intended for physicists, to provide more \"value added\" with respect to conventional programming languages."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## State of the art\n",
    "\n",
    "Two broad classes of languages have been proposed so far: ADL, CutLang, and YAML as an ADL can collectively be called languages with a \"block syntax,\" and RDataFrame, NAIL, [func-adl](https://iris-hep.org/projects/func-adl.html), and my long-defunct [Femtocode](https://github.com/diana-hep/femtocode) can be called \"Spark-like functional languages.\" Except for Femtocode, all of the above were discussed at the [Analysis Description Language Workshop](https://indico.cern.ch/event/769263/timetable/#day-2019-05-06) at Fermilab on May 6‒8, 2019.\n",
    "\n",
    "After describing the advantages and disadvantages of these languages, I'll show three toy languages that introduce new language features beyond these for various domain-specific purposes. The third of these is PartiQL."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Block languages\n",
    "\n",
    "Block languages have few recognizable programming constructs, providing a high-level feel to analysis code. There are no explicit loops, so a backend system can parallelize or distribute them however necessary. (Decisions about parallelization may still be in the user's control, but not mixed with analysis logic.) Here is an example:\n",
    "\n",
    "```\n",
    "define MR = fMR(megajets)\n",
    "define METl = met + leptonsVeto[0]\n",
    "define Rsql = sqrt(fMTR(megajets, METl) / MR)\n",
    "define Rsq = sqrt(fMTR(megajets, met) / MR)\n",
    "\n",
    "# Boost pre-selection cuts\n",
    "region preselection\n",
    "select AK4jets.size >= 3\n",
    "select AK8jets.size >= 1\n",
    "select MR > 800\n",
    "select Rsq > 0.08\n",
    "```\n",
    "\n",
    "Block constructs like `define` and `region/select` are equivalent to assignment and `if` statements in a general purpose language, but selections with meaningful physics content are named and are potentially searchable/programmatically accessible. Block languages are primarily intended to make analysis code more readable and shareable, not necessarily to make analysis code easier to write or develop. Constructs like `define`, rather than assignment, or `region/select`, rather than `if`, are not easier to use in the analysis-development process, but they can clarify the intent to another physicist.\n",
    "\n",
    "The disadvantage of block languages comes when complex processing is needed—when particles must be combined and looped over to search for solutions to various constraints (colloquially called \"combinatorics\" by physicists). Loops over particle candidates, selecting the best candidate, avoiding overlaps between candidate decays involving many particles, applying constraints like \"same flavor\" without requring a paricular flavor in the result—all of these typical analysis tasks (especially common when an analyis is in development) would significantly complicate the block structure. This kind of logic can be difficult to follow in a general programming language using `=`, `for`, and `if`, but it would get even more difficult if spread out into blocks, pushing more of the analysis description off-screen (due to length). One way out would be to allow analysis descriptions to call external code for the more complex manipulations, but too frequent use of such an \"escape valve\" would undermine the portability and preservation goals."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Functional languages\n",
    "\n",
    "Spark-like functional languages also abstract away loops over events, but in a way that doesn't give up general-purpose programming constructs. Replacing a `for` loop with a `map` or `filter` functional, with the body of the `for` loop becoming the function passed to the functional, moves the specification of which events are processed in what order from the `for` arguments into the implementation of the `map` or `filter` functional, which can be externally configured.\n",
    "\n",
    "In RDataFrame and NAIL, the functions passed to functionals are written in a general-purpose programming language: C++ (either as function pointers/references/lambda expressions or as C++ formatted strings that get wrapped as functions). In func-adl and Femtocode, on the other hand, the domain-specific language is self-embedable in the sense that functions passed into functionals are in the same language as the functionals themselves.\n",
    "\n",
    "Here is an example of passing C++ functions (as strings) to RDataFrame:\n",
    "\n",
    "```c++\n",
    "RDataFrame d(\"myTree\", \"file.root\");\n",
    "auto df = d.Define(\"p\", \"std::array<double, 4> p{px, py, pz}; return p;\")\n",
    "           .Filter(\"double p2 = 0.0; for (auto&& x : p) p2 += x*x; return sqrt(p2) < 10.0;\");\n",
    "```\n",
    "\n",
    "The advantage is that there are no restrictions on what can be computed in the loop. The disadvantage is that these functions are opaque to the processing engine: they cannot be further optimized or specialized for particle physics—they can only be executed. Functionals like this solve a single problem: the large-scale organization of the analysis. As in Vegas, what happens in the loops stays in the loops.\n",
    "\n",
    "I don't mean to minimize the importance of this one problem: it is widely recognized and important enough that frameworks built by physics collaborations have addressed it in the past. Usually, these object-oriented frameworks defined an abstract processor class with `begin`, `event`, and `end` methods for the physicist-user to override. The framework promises to call the user's `event` method on every event, just as RDataFrame calls the function passed to `Define` or `Filter` on every event. RDataFrame significantly expands the set of functionals (see [cheat sheet](https://root.cern.ch/doc/master/classROOT_1_1RDataFrame.html#cheatsheet)) well beyond `begin`, `event`, and `end`, and NAIL extends them further. However, like object-oriented frameworks, it leaves particle candidate manipulations to the general-purpose language."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "func-adl, on the other hand, lets you use the same functionals for nested collections as for collections of events, such as `Select` and `Where` in this example:\n",
    "\n",
    "```python\n",
    "jet_info = (events\n",
    "    .Select(lambda e:\n",
    "        (e.EventInfo('EventInfo'),\n",
    "         e.Jets('AntiKt4EMTopoJets'),\n",
    "         e.Tracks('InDetTrackParticles').Where(lambda t: t.pt() > 1000.0)))\n",
    "    .SelectMany(lambda e1:\n",
    "        e1[1].Select(lambda j:\n",
    "            (e1[0],j,e1[2].Where(lambda t:\n",
    "                DeltaR(t.eta(), t.phi(), j.eta(), j.phi()) < 0.2))))\n",
    "```\n",
    "\n",
    "However, these functionals, inspired by LINQ, do not make the business of particle combinations any easier than loops in a general-purpose programming language. (Femtocode, as originally conceived, wouldn't have been any better.) At best, the self-embedable feature allows for inner loops to be optimized in conjunction with the outer loops (e.g. vectorized calculations across sets of particles per event)."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Combining block languages with functional\n",
    "\n",
    "I have written three toy languages to try to influence the development of analysis description languages. The first of these, \"Jim's ADL demo,\" ([GitHub](https://nbviewer.jupyter.org/github/jpivarski/analysis-description-language/blob/master/binder/demo.ipynb), [Binder](https://mybinder.org/v2/gh/jpivarski/analysis-description-language/master?filepath=binder%2Fdemo.ipynb), October 2018) combined the block language approach with the functional to show that we can combine their strengths.\n",
    "\n",
    "This language used a curly bracket syntax to define cuts as visual blocks, but performed calculations using a functional language, like this:\n",
    "\n",
    "```\n",
    "region \"two muons\": muons.size >= 2\n",
    "{\n",
    "  zmass := muons.distincts                             # make pairs of distinct muons\n",
    "                .map(pair => (pair[0] + pair[1]))      # add the Lorentz vectors\n",
    "                .maxby(zcandidate => zcandidate.pt)    # find the maximum by pt\n",
    "                .mass                                  # but compute and return mass\n",
    "\n",
    "  count \"zmass\" by\n",
    "    regular(60, 0, 120) <- zmass\n",
    "}\n",
    "```\n",
    "\n",
    "These \"regions\" could be nested to express cut-flows and defined for multiple cuts to compute the same expression for several different cuts (a common physics task).\n",
    "\n",
    "```\n",
    "region \"slices\": true by regular(5, -5, 5) <- x\n",
    "{\n",
    "  region \"one\":   y > -3\n",
    "         \"two\":   y >  0\n",
    "         \"three\": y >  3\n",
    "  {\n",
    "    count \"counter\"\n",
    "  }\n",
    "}\n",
    "```\n",
    "\n",
    "Another basic primitive of particle physics analysis (and *not* general-purpose programming) is a systematic variation—running the same code with different input parameters. These can be blocks (nestable with `region`):\n",
    "\n",
    "```\n",
    "vary\n",
    "  \"central\":    epsilon :=  0       # assignments use := for clarity\n",
    "  \"sigma up\":   epsilon :=  0.5\n",
    "  \"sigma down\": epsilon := -0.5\n",
    "{\n",
    "  count \"histogram\" by\n",
    "    regular(100, -5, 5) <- x + epsilon\n",
    "}\n",
    "```\n",
    "\n",
    "and this nested hierarchy of blocks can be used to define a directory hierarchy of histograms:\n",
    "\n",
    "```python\n",
    "run[\"two muons\", \"zmass\"].plot()\n",
    "\n",
    "run[\"slices\", 2, \"histogram\"].plot()\n",
    "run[\"slices\", 3, \"histogram\"].plot()\n",
    "\n",
    "run[\"central\", \"histogram\"].plot()\n",
    "run[\"sigma up\", \"histogram\"].plot()\n",
    "run[\"sigma down\", \"histogram\"].plot()\n",
    "```\n",
    "\n",
    "However, the functional language implemented inside of the blocks didn't address the problem of combinatorics (any more than a general-purpose programming language does)."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Pattern matching\n",
    "\n",
    "The second toy language that I wrote demonstrated pattern matching for particle combinatorics ([code](https://github.com/diana-hep/rejig/blob/master/pattern-match/define-and-run.py), [talk](https://indico.cern.ch/event/820598/), May 2019). This language was intended to have the scope of regular expressions, but for particle physics: the user writes a tree-like description of a decay chain and associates members of predefined collections at various points in the tree. The language implementation searches for solutions that meet prescribed uniqueness and symmetry requirements.\n",
    "\n",
    "For example:\n",
    "\n",
    "```\n",
    "higgs(flavor1, flavor2) =\n",
    "    join {\n",
    "        z1 ~ {\n",
    "            lep1 ~ flavor1\n",
    "            lep2 ~ flavor1\n",
    "            mass = (lep1.p4 + lep2.p4).mass\n",
    "        }\n",
    "        z2 ~ {\n",
    "            lep1 ~ flavor2\n",
    "            lep2 ~ flavor2\n",
    "            mass = (lep1.p4 + lep2.p4).mass\n",
    "        }\n",
    "    }.filter(h => h.z1.lep1.charge != h.z1.lep2.charge and\n",
    "                  h.z2.lep1.charge != h.z2.lep2.charge)\n",
    "     .sort(h => (h.z1.mass - 91)**2 + (h.z2.mass - 91)**2)\n",
    "     [:1]\n",
    "\n",
    "higgs4e    = higgs(electrons, electrons)\n",
    "higgs4mu   = higgs(muons, muons)\n",
    "higgs2e2mu = higgs(electrons, muons)\n",
    "```\n",
    "\n",
    "defines a search for Higgs → ZZ → eeee, μμμμ, and eeμμ. The tilde operator (`~`) indicates that a variable is to be selected from a given collection (e.g. `electrons` or `muons`), which allows the user to specify flavor constraints by selecting particles from a single-flavor collection. Selections within a nested block are by default symmetric, so if `lep1` is `electron 1` and `lep2` is `electron 2`, the reverse (`lep1` is `electron 2` and `lep2` is `electron 1`) is not matched. Selections are also by default unique (across a `join`), so if `z1.lep1` is `electron 1`, then no other variable in the `join` may be `electron 1`.\n",
    "\n",
    "This language demonstrates how a syntax-bound tool can solve a common physics task, making the case that domain-specific syntax is useful for particle physics."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## PartiQL: SQL-style indexing and set operations\n",
    "\n",
    "PartiQL is the third toy language in this series (August 2019), intended to show how SQL's data model can improve data-handling in particle physics if applied to individual events, not collections of events.\n",
    "\n",
    "SQL is usually applied to large datasets, so we generally assume that if SQL or an SQL-like language is to be used in particle physics, the object of each query would be a large collection of events. I considered this several years ago and asked for suggestions from language experts [on StackOverflow](https://stackoverflow.com/q/38831961/1623645) and at StrangeLoop. It became apparent that although SQL can express physicists' queries, they work on a sufficiently specialized case that SQL only complicates the process. In particular, *every* physics query would only join data within each event, *never* across events, and this requires special `ON` clauses or tagging for `GROUP BY` to limit its attention to individual events.\n",
    "\n",
    "However, if the SQL expression is applied to each collision event in isolation, the conceptual issue is resolved. In practice, no SQL databases are optimized for this case—iterating over billions of tiny tables, rather than a few large ones—but that is a matter of implementation. Gordon Watts suggested this at the Analysis Description Language Workshop, and PartiQL is my interpretation of it.\n",
    "\n",
    "PartiQL is not SQL. There is still too much friction between the design of SQL and the practice of particle physics to demonstrate its effectiveness without many caveats. Instead, I used what I consider to be the key concepts of SQL to design a language that would be reasonable for doing physics."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Key concepts of SQL\n",
    "\n",
    "The main difference between SQL and most other programming languages is that SQL deals with sets and set operations, not lists and list iteration. Sets do not have an ordering or duplicates. (Arguably, most SQL implementations deal with \"bags,\" which do not have an ordering but may have duplicates. However, SQL without the `ALL` keyword, and the relational algebra on which it is based, is a language of sets.)\n",
    "\n",
    "Computers deal more naturally with lists because data must be serialized in some order and checking for duplicates is a non-local operation. SQL is usually implemented with an index data structure (\"keys\") that helps to maintain order-independence and uniqueness. An index also defines what equality means, particularly for `JOIN` operations that combine \"the same\" instances from different tables (with a user-supplied `ON` argument).\n",
    "\n",
    "For a (post-reconstruction) particle physics analysis, a [surrogate key index](https://en.wikipedia.org/wiki/Surrogate_key) is the most natural choice: as far as it is possible for the reconstruction software, particles are identified as distinct objects and enumerated. (Although some identified particles are fake, some real particles are undetected, and yet others are merged or double-counted, an end-user data analyst explores these possibilities with the labels received by reconstruction.)\n",
    "\n",
    "Each particle labeled by a surrogate key can then retain its identity, even if filtered, regrouped through set operations, or if new derived quantities are attached, just as SQL relations retain their identity through table manipulations, thanks to its index. This will be demonstrated with examples in later notebooks.\n",
    "\n",
    "Another key feature of SQL (in its most basic form) is that all data manipulation can be performed with basic set operations (union, intersection, difference, and their generalizations to sets of relations with `JOIN`) *instead of* general-purpose programs. SQL (up to SQL-92, but not in more recent versions of the standard) is not [Turing complete](https://increment.com/programming-languages/turing-incomplete-advantages/) and it doesn't need to be to manipulate data. Unlike the functional language we've been considering, it doesn't need functions to describe actions, and eliminating functions makes it easier to control analysis code. (Functions as substitutions, to avoid repetitiveness, can be achieved without essential complexity by adding [syntactic macros](https://en.wikipedia.org/wiki/Macro_(computer_science)#Syntactic_macros).)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Key features of PartiQL\n",
    "\n",
    "PartiQL combines the following features to provide a living example of how a domain-specific language can be more useful for particle physics than a general-purpose one.\n",
    "\n",
    "   1. Block structure for named, nested cuts and systematic variations (as in my October 2018 language).\n",
    "   2. Histograms as side-effects within the block structure (as in October 2018).\n",
    "   3. Object selection with and without replacement (as in my May 2019 language, but without the pattern-matching).\n",
    "   4. Identity defined by surrogate keys (as in SQL).\n",
    "   5. Manipulation in terms of set operations, not functions (as in SQL).\n",
    "   6. Syntactic macros to avoid repetitive typing (as in May 2019).\n",
    "\n",
    "As a demonstrator, PartiQL adheres so closely to the idea that everything can be done with set operation that no means of discovering the order of a list is provided."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Example of PartiQL\n",
    "\n",
    "Benchmark 8 from the [ADL Usability Benchmarks](https://github.com/iris-hep/adl-benchmarks-index#functionality-benchmarks) provides a good example because it is typical of physics analysis work and would require many loops and constraints in a general purpose programming language:\n",
    "\n",
    "> \"For events with at least three leptons and a same-flavor opposite-sign lepton pair, find the same-flavor opposite-sign lepton pair with the mass closest to 91.2 GeV and plot the pT of the leading other lepton.\"\n",
    "\n",
    "(I simplified the quantity to plot—numerical mathematics—because the structure manipulation is the focus of this work.)\n",
    "\n",
    "In PartiQL, the query becomes:\n",
    "\n",
    "```\n",
    "leptons = electrons union muons\n",
    "\n",
    "cut count(leptons) >= 3 named \"three_leptons\" {\n",
    "\n",
    "    Z = electrons as (lep1, lep2) union muons as (lep1, lep2)\n",
    "            where lep1.charge != lep2.charge\n",
    "            min by abs(mass(lep1, lep2) - 91.2)\n",
    "\n",
    "    third = leptons except [Z.lep1, Z.lep2] max by pt\n",
    "\n",
    "    hist third.pt by regular(100, 0, 250) named \"third_pt\" titled \"Leading other lepton pT\"\n",
    "}\n",
    "```\n",
    "\n",
    "We can freely form unions of electron and muon collections because they retain their identities, even when mixed. Pairs `(lep1, lep2)` must be drawn from either the `electrons` or the `muons` (not `leptons`) because we must ensure that `lep1` and `lep2` have the same flavor as each other, but then we can take their union because we don't care what that joint flavor is. The `(lep1, lep2, ...)` syntax performs a self-join without replacement; `electrons as lep1 cross electrons as lep2` would perform a self-join with replacement. The meaning of the `where` and `min by` clauses should be obvious.\n",
    "\n",
    "Since `Z.lep1` and `Z.lep2` retain their identity and were derived from the same collections as `leptons`, we can perform an `except` (set difference) operation to get the non-Z leptons. The `hist` keyword makes a histogram in the `cut` block that can be retrieved as\n",
    "\n",
    "```python\n",
    "histograms[\"three_leptons/third_pt\"]\n",
    "```\n",
    "\n",
    "The other notebooks in this directory go into detail about the design and implementation of PartiQL, and are executable so that you can experiment with it."
   ]
  }
 ],
 "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.7.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}