svn merge -r2:4 -r6:10 -r11:12 -r13:14 -r15:18 bhsh

ATOOLS/Math/Algebra_Interpreter.C

@@ -346,6 +346,7 @@ DEFINE_TWO_TERM_FUNCTION(Vec4D_PPerpR,"PPerpR",PPerp)
 DEFINE_TWO_TERM_FUNCTION(Vec4D_ThetaR,"ThetaR",Theta)
 DEFINE_TWO_TERM_FUNCTION(Vec4D_DEta,"DEta",DEta)
 DEFINE_TWO_TERM_FUNCTION(Vec4D_DPhi,"DPhi",DPhi)
+DEFINE_TWO_TERM_FUNCTION(Vec4D_DR,"DR",DR)
 
 Interpreter_Function::~Interpreter_Function() 
 {
@@ -422,8 +423,8 @@ DEFINE_INTERPRETER_FUNCTION(Interprete_Function)
   for (Algebra_Interpreter::Function_Map::const_reverse_iterator 
 	 fit=p_interpreter->Functions().rbegin();
        fit!=p_interpreter->Functions().rend();++fit) {
-    if ((pos=expr.rfind(fit->second->Tag()))!=std::string::npos &&
-	pos<rem) {
+    if ((pos=expr.rfind(fit->second->Tag()+'('))!=
+	std::string::npos && pos<rem) {
       func=fit->second;
       rem=pos;
     }}
@@ -648,6 +649,7 @@ Algebra_Interpreter::Algebra_Interpreter(const bool standard):
   AddFunction(new Vec4D_ThetaR());
   AddFunction(new Vec4D_DEta());
   AddFunction(new Vec4D_DPhi());
+  AddFunction(new Vec4D_DR());
 }
 
 Algebra_Interpreter::~Algebra_Interpreter()

ATOOLS/Math/Term.C

@@ -614,6 +614,7 @@ namespace ATOOLS {
   DEFINE_BINARY_VTERM_FUNCTION(Theta)
   DEFINE_BINARY_VTERM_FUNCTION(DEta)
   DEFINE_BINARY_VTERM_FUNCTION(DPhi)
+  DEFINE_BINARY_VTERM_FUNCTION(DR)
 
 }// end of namespace ATOOLS
 

ATOOLS/Math/Term.H

@@ -86,6 +86,7 @@ namespace ATOOLS {
     Term *Theta(const Term &ref) const;
     Term *DEta(const Term &ref) const;
     Term *DPhi(const Term &ref) const;
+    Term *DR(const Term &ref) const;
 
     inline void SetTag(const std::string &tag) { m_tag=tag; }
     inline void SetId(const size_t &id)        { m_id=id;   }

AddOns/Analysis/Scripts/Combine_Analysis.C

@@ -31,8 +31,19 @@ void PrintInfo()
   cout<<"          -rescale=<factor> rescale output by <factor>"<<endl;
 }
 
+#include "ATOOLS/Org/Exception.H"
+#include "ATOOLS/Org/CXXFLAGS.H"
+#ifdef USING__MPI
+#include "mpi.h"
+#endif
+
 int main(int argc,char **argv)
 {
+#ifdef USING__MPI
+  MPI::Init(argc,argv);
+#endif
+  ATOOLS::exh = new Exception_Handler();
+  ATOOLS::msg = new Message();
   if (argc<3){
     PrintInfo();
     return 0;
@@ -205,6 +216,11 @@ int main(int argc,char **argv)
   }
   if (sc>0) cout<<"successfully combined "<<inlist.size()
 		<<" directories containing "<<sc<<" histograms"<<endl; 
+  delete ATOOLS::msg;
+  delete ATOOLS::exh;
+#ifdef USING__MPI
+  MPI::Finalize();
+#endif
   return 0;
 }
  

COMIX/Main/Single_Process.C

@@ -167,6 +167,11 @@ bool COMIX::Single_Process::MapProcess()
 	p_mapproc=p_map=(*p_umprocs)[i];
 	m_oew=p_map->p_bg->MaxOrderEW();
 	m_oqcd=p_map->p_bg->MaxOrderQCD();
+	if (p_map->p_kpterms) {
+	  p_kpterms = new KP_Terms(p_map,0);
+	  p_kpterms->SetAlpha(p_map->p_bg->DInfo()->AMax());
+	  m_wgtinfo.AddMEweights(18);
+	}
 	msg_Tracking()<<"Mapped '"<<m_name<<"' -> '"<<mapname<<"'.\n";
 	std::string mapfile(rpa->gen.Variable("SHERPA_CPP_PATH")
 			    +"/Process/Comix/"+m_name+".map");
@@ -202,6 +207,11 @@ bool COMIX::Single_Process::MapProcess()
 		   <<Name()<<"' -> '"<<(*p_umprocs)[i]->Name()<<"'\n";
     if (p_bg->Map(*(*p_umprocs)[i]->p_bg,m_fmap)) {
       p_mapproc=p_map=(*p_umprocs)[i];
+      if (p_kpterms) {
+	delete p_kpterms;
+	p_kpterms = new KP_Terms(p_map,0);
+	p_kpterms->SetAlpha(p_map->p_bg->DInfo()->AMax());
+      }
       mapname=p_map->Name();
       msg_Tracking()<<"Mapped '"<<m_name<<"' -> '"
 		    <<mapname<<"'."<<std::endl;
@@ -268,7 +278,7 @@ double COMIX::Single_Process::SetZero()
     const NLO_subevtlist &rsubs(p_map?m_subs:p_bg->SubEvts());
     for (size_t i(0);i<rsubs.size();++i) rsubs[i]->Reset();
   }
-  return m_lastxs=m_last[1]=m_last[0]=0.0;
+  return m_w=m_dxs=m_lastxs=m_last[1]=m_last[0]=0.0;
 }
 
 double COMIX::Single_Process::Partonic
@@ -276,7 +286,7 @@ double COMIX::Single_Process::Partonic
 {
   Single_Process *sp(p_map!=NULL?p_map:this);
   if (mode==1 && !sp->p_scale->Scale2())
-    return m_lastxs=m_dxs+sp->GetKPTerms(m_flavs,mode);
+    return m_lastxs=m_dxs+m_w*GetKPTerms(m_flavs,mode);
   if (m_zero || !Selector()->Result()) return m_lastxs;
   for (size_t i(0);i<m_nin+m_nout;++i) {
     m_p[i]=p[i];
@@ -285,6 +295,7 @@ double COMIX::Single_Process::Partonic
   }
   if (p_map!=NULL && m_lookup && p_map->m_lookup) {
     m_dxs=p_map->m_dxs;
+    m_w=p_map->m_w;
     if (m_pinfo.m_fi.NLOType()&nlo_type::rsub) {
       const NLO_subevtlist &rsubs(p_map->p_bg->SubEvts());
       for (size_t i(0);i<rsubs.size();++i)
@@ -309,9 +320,9 @@ double COMIX::Single_Process::Partonic
       }
     }
   }
-  double kpterms(sp->GetKPTerms(m_flavs,mode));
+  double kpterms(m_w*GetKPTerms(m_flavs,mode));
   if (m_wgtinfo.m_nx) {
-    sp->FillMEWeights(m_wgtinfo);
+    FillMEWeights(m_wgtinfo);
     m_wgtinfo*=m_w;
     m_wgtinfo.m_w0=m_dxs;
   }
@@ -323,20 +334,22 @@ double COMIX::Single_Process::GetKPTerms
 {
   if (!(m_pinfo.m_fi.NLOType()&nlo_type::vsub)) return 0.0;
   if (mode==0) {
+    Single_Process *sp(p_map!=NULL?p_map:this);
     m_x[0]=m_x[1]=1.0;
     double eta0=1.0, eta1=1.0;
-    double w=p_bg->Coupling(0)/(2.0*M_PI);
+    double w=sp->p_bg->Coupling(0)/(2.0*M_PI);
+    bool map(p_map!=NULL && m_lookup && p_map->m_lookup);
     if (m_flavs[0].Strong()) {
       eta0=p_int->ISR()->X1();
-      m_x[0]=eta0+ran->Get()*(1.0-eta0);
+      m_x[0]=map?p_map->m_x[0]:eta0+ran->Get()*(1.0-eta0);
       w*=(1.0-eta0);
     }
     if (m_flavs[1].Strong()) {
       eta1=p_int->ISR()->X2();
-      m_x[1]=eta1+ran->Get()*(1.-eta1);
+      m_x[1]=map?p_map->m_x[1]:eta1+ran->Get()*(1.-eta1);
       w*=(1.0-eta1);
     }
-    p_kpterms->SetDSij(p_bg->DSij());
+    p_kpterms->SetDSij(sp->p_bg->DSij());
     p_kpterms->Calculate
       (p_int->Momenta(),m_x[0],m_x[1],eta0,eta1,w);
   }
@@ -349,14 +362,14 @@ double COMIX::Single_Process::GetKPTerms
     eta0=p_int->Momenta()[0].PPlus()/rpa->gen.PBeam(1).PMinus();
     eta1=p_int->Momenta()[1].PMinus()/rpa->gen.PBeam(0).PPlus();
   }
-  return m_w*p_kpterms->Get(m_x[0],m_x[1],eta0,eta1,fl,mode);
+  return p_kpterms->Get(m_x[0],m_x[1],eta0,eta1,fl,mode);
 }
 
 void COMIX::Single_Process::FillMEWeights(ME_wgtinfo &wgtinfo) const
 {
   wgtinfo.m_y1=m_x[0];
   wgtinfo.m_y2=m_x[1];
-  p_bg->FillMEWeights(wgtinfo);
+  (p_map?p_map:this)->p_bg->FillMEWeights(wgtinfo);
   if (p_kpterms) p_kpterms->FillMEwgts(wgtinfo);
 }
 

Manual/Commands.texi

@@ -73,6 +73,9 @@ Do not create result directory, see @ref{RESULT_DIRECTORY}.
 @item -b
 Switch to non-batch mode, see @ref{BATCH_MODE}.
 
+@item -V
+Print extended version information at startup.
+
 @item -v, --version
 Print versioning information.
 

SHERPA/Initialization/Initialization_Handler.C

@@ -926,6 +926,7 @@ void Initialization_Handler::ExtractCommandLineParameters(int argc,char * argv[]
     else if (ExtractValArg(helpsv,oit,"-l","LOG_FILE"));
     else if (ExtractValArg(helpsv,oit,"-j","PG_THREADS"));
     else if (ExtractValArg(helpsv,oit,"-g","GENERATE_RESULT_DIRECTORY","0"));
+    else if (ExtractValArg(helpsv,oit,"-V","PRINT_VERSION_INFO","1"));
     else if (par=="--version" || par=="-v"){
       msg_Out()<<"Sherpa Version "<<SHERPA_VERSION<<"."<<SHERPA_SUBVERSION<<endl;
       exit(0);
@@ -952,6 +953,7 @@ void Initialization_Handler::ExtractCommandLineParameters(int argc,char * argv[]
       msg_Out()<<"\t\t-j <threads>      set number of threads <threads>"<<endl;
       msg_Out()<<"\t\t-g                do not create result directory"<<endl;
       msg_Out()<<"\t\t-b                run in non-batch mode"<<endl;
+      msg_Out()<<"\t\t-V                print version info during runtime"<<endl;
       msg_Out()<<"\t\t-v,--version      print the version number"<<endl;
       msg_Out()<<"\t\t-h,--help         print this help message\n"<<endl;
       exit(0);

SHERPA/Main/Sherpa.C

@@ -167,7 +167,7 @@ bool Sherpa::InitializeTheRun(int argc,char * argv[])
   }
 #endif
 
-  DrawLogo();
+  DrawLogo(p_inithandler->DataReader()->GetValue("PRINT_VERSION_INFO",0));
 
   if (p_inithandler->InitializeTheFramework()) {
     if (!p_inithandler->CalculateTheHardProcesses()) return false;
@@ -283,7 +283,7 @@ bool Sherpa::SummarizeRun()
   return true; 
 }
 
-void Sherpa::DrawLogo() 
+void Sherpa::DrawLogo(const int mode) 
 { 
   msg_Info()<<"-----------------------------------------------------------------------------"<<std::endl;
   if (msg->Level()>0) msg_Out()<<"-----------    Event generation run with SHERPA started .......   -----------"<<std::endl;
@@ -339,7 +339,7 @@ void Sherpa::DrawLogo()
 	    <<"                                                                             "<<std::endl
 	    <<"-----------------------------------------------------------------------------"<<std::endl
 	    <<std::endl;
-  rpa->gen.PrintSVNVersion(msg->Info());
+  rpa->gen.PrintSVNVersion(msg->Info(),mode);
   rpa->gen.AddCitation
     (0,"The complete Sherpa package is published under \\cite{Gleisberg:2008ta}.");
 }

SHERPA/Main/Sherpa.H

@@ -17,7 +17,7 @@ namespace SHERPA {
     Event_Handler          * p_eventhandler;
     Input_Output_Handler   * p_iohandler;
 
-    void DrawLogo();  
+    void DrawLogo(const int mode);  
 
     void PrepareTerminate();
 

SHERPA/Run/init_nlo.sh

 #!/bin/bash
 
-if test $# -ne 2; then 
+if test $# -lt 2; then 
   echo "usage: $0 <sherpa exe> <input file>";
   exit 1;
 fi;
@@ -14,7 +14,7 @@ if test $tp = RS; then
 fi;
 sed -e's/}(run)/  ONLY_MAPPING_FILE 1;\n}(run)/g' < $2 > $2.$tp;
 sed -e'/NLO_QCD/ d' < $2.$tp > $2.B;
-cp -r $PWD/Process/ $nt/ 2>&1;
+test -z "$3" && cp -r $PWD/Process/ $nt/ 2>&1;
 $1 -f$2.B SHERPA_CPP_PATH=$nt;
 for i in $nt/Process/Comix/*[^\)].map; do
   if grep -q x $i; then

SHERPA/Tools/RootNtuple_Reader.C

@@ -16,6 +16,10 @@
 #include "TChain.h"
 #endif
 
+#ifdef USING__MPI
+#include "mpi.h"
+#endif
+
 using namespace SHERPA;
 using namespace PHASIC;
 using namespace ATOOLS;
@@ -102,6 +106,39 @@ RootNtuple_Reader::RootNtuple_Reader(const std::string & path,const std::string
     if (spos==std::string::npos) THROW(fatal_error,"Ivalid syntax");
     size_t i(ToType<size_t>(range.substr(0,spos)));
     size_t e(ToType<size_t>(range.substr(spos+1)));
+#ifdef USING__MPI
+    exh->MPISync();
+    int size=MPI::COMM_WORLD.Get_size();
+    int rank=MPI::COMM_WORLD.Get_rank();
+    int le=e, nact=1, values[2];
+    if (size>1) {
+      if (rank==0) {
+	msg_Info()<<"MPI Analysis {\n";
+	for (int tag=1;tag<size;++tag) if (exh->MPIStat(tag)) ++nact;
+	int inc=Max(1,(int)((e-i+1)/nact));
+	e=i+inc-1;
+	msg_Info()<<"  Rank 0 analyzes "<<basename
+		  <<"["<<i<<"-"<<e<<"].\n";
+	for (int tag=1;tag<size;++tag) {
+	  if (!exh->MPIStat(tag)) continue;
+	  values[0]=i+tag*inc;
+	  values[1]=i+(tag+1)*inc-1;
+	  if (tag==nact-1) values[1]=le;
+	  MPI::COMM_WORLD.Send(&values,2,MPI::INT,tag,tag);
+	  MPI::COMM_WORLD.Recv(&values,2,MPI::INT,MPI::ANY_SOURCE,size+tag);
+	  msg_Info()<<"  Rank "<<tag<<" analyzes "
+		    <<basename<<"["<<values[0]<<"-"<<values[1]<<"].\n";
+	}
+	msg_Info()<<"}\n";
+      }
+      else {
+	MPI::COMM_WORLD.Recv(&values,2,MPI::INT,0,rank);
+	i=values[0];
+	e=values[1];
+	MPI::COMM_WORLD.Send(&values,2,MPI::INT,0,size+rank);
+      }
+    }
+#endif
     for (;i<=e;++i) {
       std::string lfile(basename+ToString(i)+suffix);
       if (FileExists(lfile)) p_vars->p_f->Add(lfile.c_str());