Pythia DIS

new file: README.md modified: dis.sh new file: src/pythia_dis.cc

dis/src/pythia_dis.cc

-// main45.cc is a part of the PYTHIA event generator.
-// Copyright (C) 2020 Torbjorn Sjostrand.
-// PYTHIA is licenced under the GNU GPL v2 or later, see COPYING for details.
-// Please respect the MCnet Guidelines, see GUIDELINES for details.
-
-// Author: Stefan Prestel <stefan.prestel@thep.lu.se>.
-
-// Keywords: LHE file; hepmc;
-
-// This program (main45.cc) illustrates how a file with HepMC3 events can be
-// generated by Pythia8. See main44.cc for how to ouput HepMC2 events instead.
-// Note: both main44.cc and main45.cc can use the same main44.cmnd input card.
-
 #include "Pythia8/Pythia.h"
 #include "Pythia8Plugins/HepMC3.h"
 #include <unistd.h>
 
 using namespace Pythia8;
 
-//==========================================================================
-
-// Example main programm to illustrate merging.
+#include "clipp.h"
+using namespace clipp;
+using std::string;
+//______________________________________________________________________________
+
+enum class mode { none, help, list, part };
+
+struct settings {
+  double      E_electron   = 10.0;  // GeV
+  double      E_ion        = 275.0; // GeV
+  std::string outfile      = "dis.hepmc";
+  int         ion_PID      = 2212;
+  int         electron_PID = 11;
+  bool        help         = false;
+  bool        success      = false;
+  double      Q2_min       = 4.0;
+  int         N_events     = 1000;
+};
+
+template <typename T>
+void print_usage(T cli, const char* argv0)
+{
+  // used default formatting
+  std::cout << "Usage:\n" << usage_lines(cli, argv0) << "\nOptions:\n" << documentation(cli) << '\n';
+}
+//______________________________________________________________________________
+
+template <typename T>
+void print_man_page(T cli, const char* argv0)
+{
+  // all formatting options (with their default values)
+  auto fmt = clipp::doc_formatting{}
+                 .start_column(8)                    // column where usage lines and documentation starts
+                 .doc_column(20)                     // parameter docstring start col
+                 .indent_size(4)                     // indent of documentation lines for children of a
+                                                     // documented group
+                 .line_spacing(0)                    // number of empty lines after single documentation lines
+                 .paragraph_spacing(1)               // number of empty lines before and after paragraphs
+                 .flag_separator(", ")               // between flags of the same parameter
+                 .param_separator(" ")               // between parameters
+                 .group_separator(" ")               // between groups (in usage)
+                 .alternative_param_separator("|")   // between alternative flags
+                 .alternative_group_separator(" | ") // between alternative groups
+                 .surround_group("(", ")")           // surround groups with these
+                 .surround_alternatives("(", ")")    // surround group of alternatives with these
+                 .surround_alternative_flags("", "") // surround alternative flags with these
+                 .surround_joinable("(", ")")        // surround group of joinable flags with these
+                 .surround_optional("[", "]")        // surround optional parameters with these
+                 .surround_repeat("", "...");        // surround repeatable parameters with these
+
+  auto mp = make_man_page(cli, argv0, fmt);
+  mp.prepend_section("DESCRIPTION", "Simple pythia dis generator.");
+  mp.append_section("EXAMPLES", " $ pythia_dis [output_file]");
+  std::cout << mp << "\n";
+}
+//______________________________________________________________________________
+
+settings cmdline_settings(int argc, char* argv[]) {
+  settings s;
+  auto     lastOpt =
+      " options:" % (option("-h", "--help").set(s.selected, mode::help) % "show help",
+                     option("-o", "--output") & value("out", s.outfile), value("file", s.infile).if_missing([] {
+                       std::cout << "You need to provide an input xml filename as the last "
+                                    "argument!\n";
+                     }) % "input xml file");
+
+  auto server_cli =
+      ((option("-p", "--port") & value("http_port", s.http_port)) %
+           "port to which the http serve attaches. Default: 8090 ",
+       (option("-H", "--host") & value("http_host", s.http_host)) %
+           "Http server host name or IP address. Default: 127.0.0.1",
+       (option("-f", "--file") & value("io_file", s.in_out_file)) %
+           "File used to initialize and save plots. Default: top_folder.root");
+
+  auto cli =
+      (command("help").set(s.selected, mode::help) | (server_cli, lastOpt));
+
+  assert(cli.flags_are_prefix_free());
+
+  auto res = parse(argc, argv, cli);
+
+  // if( res.any_error() ) {
+  //  s.success = false;
+  //  std::cout << make_man_page(cli, argv[0]).prepend_section("error: ",
+  //                                                           "    The best
+  //                                                           thing since
+  //                                                           sliced bread.");
+  //  return s;
+  //}
+  s.success = true;
+
+  if (s.selected == mode::help) {
+    print_man_page<decltype(cli)>(cli, argv[0]);
+  };
+  return s;
+}
+//______________________________________________________________________________
 
-int main( int argc, char* argv[] ){
+int main(int argc, char* argv[]) {
 
-  // Check that correct number of command-line arguments
-  if (argc != 2) {
-    cerr << " Unexpected number of command-line arguments ("<<argc<<"). \n"
-         << " You are expected to provide the arguments" << endl
-         << " 1. Output file for HepMC events" << endl
-         << " Program stopped. " << endl;
+  settings s = cmdline_settings(argc, argv);
+  if (!s.success) {
     return 1;
   }
 
-  // Beam energies, minimal Q2, number of events to generate.
-  double eProton   = 250.;
-  double eElectron = 10.0;
-  double Q2min     = 5.;
-  int    nEvent    = 10000;
+  if (s.selected == mode::help) {
+    return 0;
+  }
+
+
+    // Beam energies, minimal Q2, number of events to generate.
+  double eProton   = s.E_ion;
+  double eElectron = s.E_electron;
+  double Q2min     = s.Q2_min;
+  int    nEvent    = s.N_events;
 
   // Generator. Shorthand for event.
   Pythia pythia;
@@ -45,7 +127,7 @@ int main( int argc, char* argv[] ){
   // Set up incoming beams, for frame with unequal beam energies.
   pythia.readString("Beams:frameType = 2");
   // BeamA = proton.
-  pythia.readString("Beams:idA = 2212");
+  pythia.readString("Beams:idA = " +std::to_string(s.ion_PID));
   pythia.settings.parm("Beams:eA", eProton);
   // BeamB = electron.
   pythia.readString("Beams:idB = 11");
@@ -55,7 +137,7 @@ int main( int argc, char* argv[] ){
   // Neutral current (with gamma/Z interference).
   pythia.readString("WeakBosonExchange:ff2ff(t:gmZ) = on");
   // Uncomment to allow charged current.
-  //pythia.readString("WeakBosonExchange:ff2ff(t:W) = on");
+  // pythia.readString("WeakBosonExchange:ff2ff(t:W) = on");
   // Phase-space cut: minimal Q2 of process.
   pythia.settings.parm("PhaseSpace:Q2Min", Q2min);
 
@@ -80,21 +162,21 @@ int main( int argc, char* argv[] ){
   cout << endl << endl << endl;
 
   // Histograms.
-  double Wmax = sqrt(4.* eProton * eElectron);
-  Hist Qhist("Q [GeV]", 100, 0., 50.);
-  Hist Whist("W [GeV]", 100, 0., Wmax);
-  Hist xhist("x", 100, 0., 1.);
-  Hist yhist("y", 100, 0., 1.);
-  Hist pTehist("pT of scattered electron [GeV]", 100, 0., 50.);
-  Hist pTrhist("pT of radiated parton [GeV]", 100, 0., 50.);
-  Hist pTdhist("ratio pT_parton/pT_electron", 100, 0., 5.);
+  double Wmax = sqrt(4. * eProton * eElectron);
+  Hist   Qhist("Q [GeV]", 100, 0., 50.);
+  Hist   Whist("W [GeV]", 100, 0., Wmax);
+  Hist   xhist("x", 100, 0., 1.);
+  Hist   yhist("y", 100, 0., 1.);
+  Hist   pTehist("pT of scattered electron [GeV]", 100, 0., 50.);
+  Hist   pTrhist("pT of radiated parton [GeV]", 100, 0., 50.);
+  Hist   pTdhist("ratio pT_parton/pT_electron", 100, 0., 5.);
 
   double sigmaTotal(0.), errorTotal(0.);
-  bool wroteRunInfo = false;
-    // Get the inclusive x-section by summing over all process x-sections.
-    double xs = 0.;
-    for (int i=0; i < pythia.info.nProcessesLHEF(); ++i)
-      xs += pythia.info.sigmaLHEF(i);
+  bool   wroteRunInfo = false;
+  // Get the inclusive x-section by summing over all process x-sections.
+  double xs = 0.;
+  for (int i = 0; i < pythia.info.nProcessesLHEF(); ++i)
+    xs += pythia.info.sigmaLHEF(i);
 
   // Begin event loop.
   for (int iEvent = 0; iEvent < nEvent; ++iEvent) {
@@ -176,8 +258,6 @@ int main( int argc, char* argv[] ){
     // Write the HepMC event to file. Done with it.
     ascii_io.write_event(hepmcevt);
 
-
-
     // End of event loop. Statistics and histograms.
   }
   pythia.stat();
@@ -185,5 +265,4 @@ int main( int argc, char* argv[] ){
 
 
   return 0;
-
-}
+  }