|
||||
File indexing completed on 2025-01-18 09:58:46
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 // Author: Dennis Wright (SLAC) 0027 // Date: 28 January 2013 0028 // 0029 // Description: implementation base for numerically integrated two-body 0030 // final state angular distributions in Bertini-style cascade 0031 // 0032 // 20130219 M. Kelsey: Move to .icc file, required for templated classes 0033 // 20130620 Fix Coverity #50295, recursive #include. 0034 // 20130924 M. Kelsey -- Use G4Log, G4Exp for CPU speedup 0035 0036 #ifndef G4NumIntTwoBodyAngDst_icc 0037 #define G4NumIntTwoBodyAngDst_icc 1 0038 0039 #include "G4Log.hh" 0040 #include "G4Exp.hh" 0041 #include "Randomize.hh" 0042 0043 0044 template <G4int NKEBINS, G4int NANGLES> G4double 0045 G4NumIntTwoBodyAngDst<NKEBINS,NANGLES>:: 0046 GetCosTheta(const G4double& ekin, const G4double& pcm) const 0047 { 0048 G4double costh=1.0, slope=1.0; // Buffers for final result 0049 G4double xrand = G4UniformRand(); // Random value to sample CDF 0050 0051 if (ekin < labKE[nDists-1]) { 0052 Interpolate(ekin); // Get intermediate CDF @ energy 0053 for (G4int i=1; i<nAngles; i++) { 0054 if (xrand < angDist[i]) { 0055 slope = (cosBins[i] - cosBins[i-1])/(angDist[i] - angDist[i-1]); 0056 costh = (xrand - angDist[i-1])*slope + cosBins[i-1]; 0057 break; 0058 } 0059 } 0060 } else { 0061 slope = 2.*tcoeff*pcm*pcm; 0062 costh = G4Log(1.0 - xrand*(1. - G4Exp(2.*slope) ) )/slope - 1.0; 0063 } 0064 return costh; 0065 } 0066 0067 0068 // Find normalized indefinite integral of angular distribution at interpolated 0069 // kinetic energy 0070 0071 template <G4int NKEBINS, G4int NANGLES> 0072 void G4NumIntTwoBodyAngDst<NKEBINS,NANGLES>::Interpolate(const G4double& ekin) const 0073 { 0074 for (G4int j=1; j<nDists; j++) { 0075 if (ekin >= labKE[j]) continue; // Find desired pair of cos(theta) CDFs 0076 0077 G4double frac = (ekin - labKE[j-1])/(labKE[j] - labKE[j-1]); 0078 for (G4int i=0; i < nAngles; i++) { 0079 angDist[i] = (1.0 - frac)*angDists[j-1][i] + frac*angDists[j][i]; 0080 } 0081 0082 break; 0083 } 0084 } 0085 0086 #endif /* G4NumIntTwoBodyAngDst_icc */
[ Source navigation ] | [ Diff markup ] | [ Identifier search ] | [ general search ] |
This page was automatically generated by the 2.3.7 LXR engine. The LXR team |