From a79e267a23ad29d24550e14ca97c41a0af95c373 Mon Sep 17 00:00:00 2001 From: Felix Riehn <felix@matilda> Date: Thu, 20 May 2021 17:10:02 +0100 Subject: [PATCH] clean up --- corsika/detail/modules/epos/Interaction.inl | 59 -------------- corsika/modules/epos/Interaction.hpp | 10 +-- modules/epos/epos.hpp | 89 +-------------------- tests/modules/testEpos.cpp | 11 --- 4 files changed, 4 insertions(+), 165 deletions(-) diff --git a/corsika/detail/modules/epos/Interaction.inl b/corsika/detail/modules/epos/Interaction.inl index 31c0886c5..fc29c25be 100644 --- a/corsika/detail/modules/epos/Interaction.inl +++ b/corsika/detail/modules/epos/Interaction.inl @@ -87,11 +87,9 @@ namespace corsika::epos { // corsika7 ini int iarg = 0; ::epos::aaset_(iarg); - //::epos::atitle_(); ::epos::prnt1_.ish = 0; // debug level in epos ::epos::files_.ifch = 6; // output unit - //::epos::prnt1_.iecho = 1; // dummy set seeds for random number generator in epos. need to fool epos checks... // we will use external generator @@ -148,22 +146,6 @@ namespace corsika::epos { strcpy(::epos::fname_.fncs, CS.data); ::epos::nfname_.nfncs = CS.length; - // dummy event - // ::epos::hadr25_.idprojin = 1120; - // ::epos::hadr25_.idtargin = 1120; - // //#if __CONEX__ && __EPOS__ && __HIGHMEM__ - // // maproj = 250 - // //#else - // ::epos::nucl1_.maproj = 56; - // // #endif - // ::epos::nucl1_.laproj = 28; - // ::epos::nucl1_.matarg = 14; - // ::epos::nucl1_.latarg = 1; - // ::epos::hadr1_.pnll = 200.; - // ::epos::lept1_.engy = -1.; - - //::epos::ainit_(); - // dummy event (prepare commons) initializeEventLab(Code::Iron, Iron::nucleus_A, Iron::nucleus_Z, Code::Argon, Argon::nucleus_A, Argon::nucleus_Z, 100_GeV); @@ -285,18 +267,6 @@ namespace corsika::epos { inline Interaction::~Interaction() { CORSIKA_LOGGER_DEBUG(logger_, "n={} ", count_); } - // inline std::tuple<corsika::CrossSectionType, corsika::CrossSectionType> - // Interaction::getCrossSection(corsika::Code const BeamId, corsika::Code const - // TargetId, - // const corsika::HEPEnergyType EnergyCOM) const { - // if (!is_nucleus(BeamId)) - // return getCrossSection(BeamId, 1, 1, TargetId, get_nucleus_A(TargetId), - // get_nucleus_Z(TargetId), EnergyCOM); - // else - // throw("nuclear projecile, call getCrossSection with : BeamId, BeamA, BeamZ, - // ..."); - // } - inline std::tuple<corsika::CrossSectionType, corsika::CrossSectionType> Interaction::calcCrossSectionCoM(corsika::Code const BeamId, int const BeamA, int const BeamZ, corsika::Code const TargetId, @@ -336,8 +306,6 @@ namespace corsika::epos { "calcCrossSectionCoM: interaction of beam hadron not defined in " "Epos!"); - //::epos::xsigma_(); - double sigProd, sigEla = 0; float sigTot1, sigProd1, sigEla1, sigCut1 = 0; if (!is_nucleus(TargetId) && !is_nucleus(BeamId)) { @@ -345,16 +313,9 @@ namespace corsika::epos { sigEla = ::epos::hadr5_.sigela; } else { ::epos::crseaaepos_(sigTot1, sigProd1, sigCut1, sigEla1); - // sigProd = ::epos::hadr5_.sigineaa; - // sigEla = ::epos::hadr5_.sigelaaa; sigProd = sigProd1; sigEla = sigEla1; } - - // calculate cross section - // float sigTot, sigProd, sigEla, sigCut = 0; - //::epos::crseaaepos_(sigTot, sigProd, sigCut, sigEla); - CORSIKA_LOGGER_DEBUG(logger_, "calcCrossSectionCoM: output:" " sigProd={} mb," @@ -465,12 +426,6 @@ namespace corsika::epos { } // get target from environment - /* - the target should be defined by the Environment, - ideally as full particle object so that the four momenta - and the boosts can be defined.. - */ - MomentumVector const& pLab = projectile.getMomentum(); CoordinateSystemPtr const& labCS = pLab.getCoordinateSystem(); @@ -479,14 +434,6 @@ namespace corsika::epos { // total momentum and energy HEPEnergyType Elab = projectile.getEnergy() + constants::nucleonMass; - // MomentumVector pTotLab(labCS, {0_GeV, 0_GeV, 0_GeV}); - // pTotLab += pLab; - // pTotLab += pTarget; - // auto const pTotLabNorm = pTotLab.getNorm(); - // calculate cm. energy - // const HEPEnergyType ECoM = sqrt( - // (Elab + pTotLabNorm) * (Elab - pTotLabNorm)); // binomial for numerical - // accuracy auto const* currentNode = projectile.getNode(); const auto& mediumComposition = @@ -548,11 +495,6 @@ namespace corsika::epos { HEPEnergyType const projectileMomentumLabPerNucleon = projectileMomentum / beamA; // define target - // auto const eTargetLab = 0_GeV + constants::nucleonMass; - // HEPEnergyType const Etot = eProjectileLab + eTargetLab; - // invariant mass, i.e. cm. energy - // HEPEnergyType const Ecm = - // sqrt((Etot + projectileMomentum) * (Etot - projectileMomentum)); // sample target mass number auto const* currentNode = projectile.getNode(); @@ -581,7 +523,6 @@ namespace corsika::epos { // // from corsika7 interface // // NEXLNK-part - int targetA = 1; int targetZ = 1; if (is_nucleus(targetCode)) { diff --git a/corsika/modules/epos/Interaction.hpp b/corsika/modules/epos/Interaction.hpp index ec6cbeeaf..2288fdebd 100644 --- a/corsika/modules/epos/Interaction.hpp +++ b/corsika/modules/epos/Interaction.hpp @@ -45,19 +45,12 @@ namespace corsika::epos { Code const, int const, int const, Code const, int const, int const, HEPEnergyType const) const; - // std::tuple<CrossSectionType, CrossSectionType> getCrossSection( - // Code const, int const, int const, Code const, int const, int const, - // HEPEnergyType const) const; - - // std::tuple<CrossSectionType, CrossSectionType> getCrossSection( - // Code const, Code const, HEPEnergyType const) const; - template <typename TParticle> GrammageType getInteractionLength(TParticle const&) const; /** In this function EPOSLHC is called to produce one event. The - event is copied (and boosted) into the shower lab frame. + event is copied into the shower lab frame. */ template <typename TSecondaries> void doInteraction(TSecondaries&); @@ -65,6 +58,7 @@ namespace corsika::epos { bool isValidCoMEnergy(HEPEnergyType const ecm) const { return (minEnergyCoM_ <= ecm) && (ecm <= maxEnergyCoM_); } + //! eposlhc only accepts nuclei with X<=A<=Y as targets, or protons aka Hydrogen or //! neutrons (p,n == nucleon) bool isValidTarget(Code const) const; diff --git a/modules/epos/epos.hpp b/modules/epos/epos.hpp index 7419bf1b2..48809a608 100644 --- a/modules/epos/epos.hpp +++ b/modules/epos/epos.hpp @@ -52,8 +52,6 @@ namespace epos { void conini_(); void psaini_(); - // void idspin_(int&, int&, int&, int&); - // void iclass_(int&, int&); void emsini_(double&, int&, int&); void paramini_(int&); void xsigma_(); @@ -87,8 +85,6 @@ namespace epos { // c------------------------------------------------------------------------------ void crseaaepos_(float&, float&, float&, float&); - // double cxepocrse_(double&, int&, int&, int&); - void emsaaa_(int&); void gakfra_(int&, int&); void utghost_(int&); @@ -108,13 +104,10 @@ namespace epos { // convert id from one format to another int idtrafo_(char[3], char[3], int&); - // conex routine - // void cxidmass_(int&, int&); - // additional random number functions void ranfini_(double&, int&, int&); void ranfcv_(double&); - // void ranfgt(int&); + float rangen_(); double drangen_(); // common blocks as @@ -185,19 +178,13 @@ namespace epos { int jdecay; int iremn; } othe2_; - // integer ifrade,iframe,idecay,jdecay,iremn - // common/othe2/ifrade,iframe,idecay,jdecay,iremn extern struct { int ktnbod; } metr7_; - // integer ktnbod - // common/metr7/ktnbod extern struct { float egylow; float egyfac; } had12_; - // real egylow,egyfac - // common/had12/egylow,egyfac extern struct { int laproj; @@ -207,9 +194,6 @@ namespace epos { float core; float fctrmx; } nucl1_; - // real core,fctrmx - // integer laproj,maproj,latarg,matarg - // common/nucl1/laproj,maproj,latarg,matarg,core,fctrmx extern struct { float amproj; @@ -218,8 +202,6 @@ namespace epos { float yhaha; float pnullx; } chadron_; - // real amproj,amtarg,ypjtl,yhaha,pnullx - // common/chadron/amproj,amtarg,ypjtl,yhaha,pnullx extern struct { int iomodl; @@ -227,9 +209,6 @@ namespace epos { int idtarg; float wexcit; } hadr2_; - // integer iomodl,idproj,idtarg - // real wexcit - // common/hadr2/iomodl,idproj,idtarg,wexcit extern struct { int idprojin; @@ -241,11 +220,6 @@ namespace epos { int isotarg; } hadr25_; - // real rexdifi,rexndii - // integer idprojin,idtargin,irdmpr,isoproj,isotarg - // common/hadr25/idprojin,idtargin,rexdifi(4),rexndii(4),irdmpr, - // * isoproj,isotarg - extern struct { float engy; float elepti; @@ -253,10 +227,6 @@ namespace epos { float angmue; int icinpu; } lept1_; - // real engy, elepti, elepto, angmue integer icinpu - // common / lept1 / engy, elepti, elepto, - // angmue, - // icinpu extern struct { float egymin; @@ -265,9 +235,6 @@ namespace epos { float ecms; float ekin; } enrgy_; - // real egymin, egymax, elab, ecms, ekin - // common / enrgy / egymin, egymax, elab, ecms, - // ekin extern struct { float pnll; @@ -277,10 +244,6 @@ namespace epos { float wproj; float wtarg; } hadr1_; - // real pnll, ptq, exmass, cutmss, wproj, wtarg - // common / hadr1 / pnll, ptq, exmass, cutmss, - // wproj, - // wtarg unsigned int constexpr idxD0 = 0; unsigned int constexpr idxD1 = 2; @@ -301,17 +264,6 @@ namespace epos { float bmxdif[nclha][nclha]; float bkmxndif; } Dparam_; - // real bmxdif,bkmxndif - // integer idxDmin - // common / Dparam / alpD(idxD0: idxD1, nclha, nclha), - //* alpDp(idxD0 : idxD1, nclha, nclha), - //*alpDpp(idxD0 : idxD1, nclha, nclha), - //* betD(idxD0 : idxD1, nclha, nclha), - //* betDp(idxD0 : idxD1, nclha, nclha), - //*betDpp(idxD0 : idxD1, nclha, nclha), - //* gamD(idxD0 : idxD1, nclha, nclha), - //* delD(idxD0 : idxD1, nclha, nclha), - //*idxDmin, bmxdif(nclha, nclha), bkmxndif extern struct { float phievt; @@ -340,16 +292,7 @@ namespace epos { int maxfra; int kohevt; } cevt_; - // real phievt, bimevt, pmxevt, egyevt , xbjevt, qsqevt, zppevt, zptevt - // integer nevt, - // kolevt, koievt, kohevt, npjevt , ntgevt, npnevt, nppevt, ntnevt, ntpevt, jpnevt, - // jppevt, jtnevt, jtpevt , nglevt, minfra, maxfra - // common / cevt / phievt, nevt, - // bimevt, kolevt, koievt, pmxevt, egyevt, npjevt , ntgevt, npnevt, nppevt, ntnevt, - // ntpevt, jpnevt, jppevt, jtnevt, jtpevt , xbjevt, qsqevt, nglevt, zppevt, zptevt, - // minfra, maxfra, kohevt - - // common/cseed/seedi,seedj,seedj2,seedc,iseqini,iseqsim + extern struct { double seedi; double seedj; @@ -359,8 +302,6 @@ namespace epos { int iseqsim; } cseed_; - // integer istore,istmax,irescl,ntrymx,nclean,iopdg,ioidch - // common/othe1/istore,istmax,gaumx,irescl,ntrymx,nclean,iopdg,ioidch extern struct { int istore; int istmax; @@ -372,8 +313,6 @@ namespace epos { int ioidch; } othe1_; - // integer ifop,ifmt,ifch,ifcx,ifhi,ifdt,ifcp,ifdr - // common/files/ifop,ifmt,ifch,ifcx,ifhi,ifdt,ifcp,ifdr extern struct { int ifop; int ifmt; @@ -385,10 +324,6 @@ namespace epos { int ifdr; } files_; - // character*500 fnch,fnhi,fndt,fnii,fnid,fnie,fnrj,fnmt - // * ,fngrv,fncp,fnnx,fncs,fndr,fnhpf - // common/fname/ fnch, fnhi, fndt, fnii, fnid, fnie, fnrj, fnmt - // * ,fngrv,fncp,fnnx,fncs,fndr,fnhpf extern struct { char fnch[500]; char fnhi[500]; @@ -407,10 +342,6 @@ namespace epos { } fname_; - // integer nfnch,nfnhi,nfndt,nfnii,nfnid,nfnie,nfnrj,nfnmt - // *,nfngrv,nfncp,nfnnx,nfncs,nfndr,nfnhpf - // common/nfname/nfnch,nfnhi,nfndt,nfnii,nfnid,nfnie,nfnrj,nfnmt - // *,nfngrv,nfncp,nfnnx,nfncs,nfndr,nfnhpf extern struct { int nfnch; int nfnhi; @@ -428,8 +359,6 @@ namespace epos { int nfnhpf; } nfname_; - // integer iprmpt,ish,ishsub,irandm,irewch,iecho,modsho,idensi - // common/prnt1/iprmpt,ish,ishsub,irandm,irewch,iecho,modsho,idensi extern struct { int iprmpt; int ish; @@ -442,11 +371,6 @@ namespace epos { } prnt1_; unsigned int constexpr mmry = 1; unsigned int constexpr mxptl = 200000 / mmry; - // real pptl,tivptl,xorptl - // integer nptl,iorptl,idptl,istptl,ifrptl,jorptl,ibptl,ityptl - // common/cptl/nptl,pptl(5,mxptl),iorptl(mxptl),idptl(mxptl) - // *,istptl(mxptl),tivptl(2,mxptl),ifrptl(2,mxptl),jorptl(mxptl) - // *,xorptl(4,mxptl),ibptl(4,mxptl),ityptl(mxptl) extern struct { int nptl; float pptl[mxptl][5]; @@ -461,10 +385,6 @@ namespace epos { int ityptl[mxptl]; } cptl_; - // real sigtot,sigcut,sigela,sloela,sigsd,sigine,sigdif - // *,sigineaa,sigtotaa,sigelaaa,sigcutaa,sigdd - // common/hadr5/sigtot,sigcut,sigela,sloela,sigsd,sigine,sigdif - // *,sigineaa,sigtotaa,sigelaaa,sigcutaa,sigdd extern struct { float sigtot; float sigcut; @@ -480,17 +400,12 @@ namespace epos { float sigdd; } hadr5_; - // integer mxnody,nrnody,nody - // parameter(mxnody=200) - // common/nodcy/nrnody,nody(mxnody) unsigned int constexpr mxnody = 200; extern struct { int nrnody; int nody[mxnody]; } nodcy_; - // integer iclpro,icltar,iclegy - // common/had10/iclpro,icltar,iclegy extern struct { int iclpro; int icltar; diff --git a/tests/modules/testEpos.cpp b/tests/modules/testEpos.cpp index c5ed8743c..282b9d178 100644 --- a/tests/modules/testEpos.cpp +++ b/tests/modules/testEpos.cpp @@ -89,12 +89,6 @@ TEST_CASE("Epos", "[processes]") { } } } - - // SECTION("validation - mass") { - // for (auto p : get_all_particles()) - // if (!is_nucleus(p)) - // CHECK(get_mass(p) / corsika::epos::getEposMass(p) == Approx(1).margin(0.5)); - // } } #include <corsika/framework/geometry/Point.hpp> @@ -144,11 +138,8 @@ TEST_CASE("EposInterface", "[processes]") { SECTION("InteractionInterface - valid targets") { Interaction model; - // eposlhc accepts protons or nuclei with 4<=A<=18 as targets CHECK_FALSE(model.isValidTarget(Code::Electron)); CHECK(model.isValidTarget(Code::Hydrogen)); - // CHECK_FALSE(model.isValidTarget(Code::Deuterium)); - // CHECK_FALSE(model.isValidTarget(Code::Helium3)); CHECK(model.isValidTarget(Code::Helium)); CHECK_FALSE(model.isValidTarget(Code::Iron)); CHECK(model.isValidTarget(Code::Oxygen)); @@ -169,7 +160,6 @@ TEST_CASE("EposInterface", "[processes]") { SECTION("InteractionInterface - hadron cross sections") { Interaction model; - // eposlhc accepts protons or nuclei with 4<=A<=18 as targets // p-p at 7TeV around 70mb according to LHC auto const [xs_prod, xs_ela] = @@ -193,7 +183,6 @@ TEST_CASE("EposInterface", "[processes]") { SECTION("InteractionInterface - nuclear cross sections") { Interaction model; - // eposlhc accepts protons or nuclei with 4<=A<=18 as targets auto const [xs_prod, xs_ela] = model.getCrossSectionLab( Code::Proton, 1, 1, Code::Oxygen, Oxygen::nucleus_A, Oxygen::nucleus_Z, 100_GeV); -- GitLab