include/Geant4/G4CrossSectionBuffer.hh

0001 //
0002 // ********************************************************************
0003 // * License and Disclaimer                                           *
0004 // *                                                                  *
0005 // * The  Geant4 software  is  copyright of the Copyright Holders  of *
0006 // * the Geant4 Collaboration.  It is provided  under  the terms  and *
0007 // * conditions of the Geant4 Software License,  included in the file *
0008 // * LICENSE and available at  http://cern.ch/geant4/license .  These *
0009 // * include a list of copyright holders.                             *
0010 // *                                                                  *
0011 // * Neither the authors of this software system, nor their employing *
0012 // * institutes,nor the agencies providing financial support for this *
0013 // * work  make  any representation or  warranty, express or implied, *
0014 // * regarding  this  software system or assume any liability for its *
0015 // * use.  Please see the license in the file  LICENSE  and URL above *
0016 // * for the full disclaimer and the limitation of liability.         *
0017 // *                                                                  *
0018 // * This  code  implementation is the result of  the  scientific and *
0019 // * technical work of the GEANT4 collaboration.                      *
0020 // * By using,  copying,  modifying or  distributing the software (or *
0021 // * any work based  on the software)  you  agree  to acknowledge its *
0022 // * use  in  resulting  scientific  publications,  and indicate your *
0023 // * acceptance of all terms of the Geant4 Software license.          *
0024 // ********************************************************************
0025 //
0026 #ifndef G4CrossSectionBuffer_h
0027 #define G4CrossSectionBuffer_h
0028
0029 #include <utility>
0030 #include <vector>
0031 #include <CLHEP/Units/SystemOfUnits.h>
0032
0033 #include "globals.hh"
0034 #include "G4ParticleDefinition.hh"
0035 #include "G4KineticTrack.hh"
0036
0037 class G4CrossSectionBuffer
0038 {
0039   public:
0040
0041     G4CrossSectionBuffer(const G4ParticleDefinition * aA, const G4ParticleDefinition * aB)
0042     : theA(aA), theB(aB) {}
0043
0044     G4bool InCharge(const G4ParticleDefinition * aA, const G4ParticleDefinition * aB) const
0045     {
0046       G4bool result = false;
0047       if(aA == theA)
0048       {
0049         if(aB == theB) result = true;
0050       }
0051       else if(aA == theB)
0052       {
0053         if(aB == theA) result = true;
0054       }
0055       return result;
0056     }
0057
0058     void push_back(G4double S, G4double x)
0059     {
0060       std::pair<G4double, G4double> aNew;
0061       aNew.first = S;
0062       aNew.second = x;
0063       theData.push_back(aNew);
0064     }
0065     G4double CrossSection(const G4KineticTrack& trk1, const G4KineticTrack& trk2) const
0066     {
0067       G4double sqrts = (trk1.Get4Momentum()+trk2.Get4Momentum()).mag();
0068       G4double x1(1), y1(0);
0069       G4double x2(2), y2(0);
0070
0071       if(theData.size()==1) return theData[theData.size()-1].second;
0072
0073       for(size_t i=0; i<theData.size(); i++)
0074       {
0075         if(theData[i].first>sqrts)
0076     {
0077       if(0==i)
0078       {
0079         x1 = theData[i].first;
0080         y1 = theData[i].second;
0081         x2 = theData[i+1].first;
0082         y2 = theData[i+1].second;
0083       }
0084       else if(theData.size()-1==i)
0085       {
0086         x1 = theData[theData.size()-2].first;
0087         y1 = theData[theData.size()-2].second;
0088         x2 = theData[theData.size()-1].first;
0089         y2 = theData[theData.size()-1].second;
0090       }
0091       else
0092       {
0093         x1 = theData[i-1].first;
0094         y1 = theData[i-1].second;
0095         x2 = theData[i].first;
0096         y2 = theData[i].second;
0097       }
0098       break;
0099     }
0100       }
0101       // LINLIN interpolation
0102       // G4cerr << "!!!!!! "<<sqrts<<" "<<x1<<" "<<x2<<" "<<y1<<" "<<y2<<" "<<i<<G4endl;
0103       G4double result = y1 + (sqrts-x1) * (y2-y1)/(x2-x1);
0104       if(result<0) result = 0;
0105       if(y1<0.01*CLHEP::millibarn) result = 0;
0106       return result;
0107     }
0108
0109     void Print()
0110     {
0111       for(size_t i=0;i<theData.size(); i++)
0112       {
0113         G4cerr << "sqrts = "<<theData[i].first<<", X = "<<theData[i].second/CLHEP::millibarn<<G4endl;
0114       }
0115     }
0116
0117   private:
0118   std::vector<std::pair<G4double, G4double> > theData;
0119
0120   const G4ParticleDefinition * theA;
0121   const G4ParticleDefinition * theB;
0122 };
0123
0124 #endif