Skip to content
Snippets Groups Projects
Commit 5529b220 authored by Bertram Kopf's avatar Bertram Kopf
Browse files

added dynamics for sigma parametrization

parent a95d18a3
No related branches found
No related tags found
No related merge requests found
......@@ -73,7 +73,10 @@ complex<double> SigmaParameterization::calc(double currentMass){
return result;
}
std::complex<double> SigmaParameterization::php4pi(double currentMass2){
complex<double> SigmaParameterization::calc(double currentMass, double gSigma){
return gSigma*calc(currentMass);
}
std::complex<double> SigmaParameterization::php4pi(double currentMass2){
// if(currentMass2<16.*PawianConstants::mPiSq){
// Alert << "current mass^2 must be larger than 16*m_pi^2!!!\n current mass^2: " << currentMass2 <<"\n16*m_pi^2: " << 16.*PawianConstants::mPiSq << endmsg;
// exit(1);
......
......@@ -45,6 +45,7 @@ public:
~SigmaParameterization();
complex<double> calc(double currentMass);
complex<double> calc(double currentMass, double gSigma);
complex<double> calcT(double currentMass);
protected:
......
......@@ -56,6 +56,8 @@
#include "PwaUtils/LinearDynamics.hh"
#include "PwaUtils/OmnesDynamics.hh"
#include "PwaUtils/SExpDynamics.hh"
#include "PwaUtils/SigmaParamDynamics.hh"
#include "PwaUtils/ProdChannelInfo.hh"
#include "PwaUtils/GlobalEnv.hh"
......@@ -190,6 +192,9 @@ std::shared_ptr<AbsDynamics> DynRegistry::getDynamics(std::shared_ptr<AbsDecay>
}
result= std::shared_ptr<AbsDynamics>(new BreitWignerBlattWTensorRelDynamics(theName, fsParticles, theDec->motherPart(), fsParticlesDaughter1, fsParticlesDaughter2, theDec->barrierqR()));
}
else if(theDec->dynType()=="SigmaParam"){
result= std::shared_ptr<AbsDynamics>(new SigmaParamDynamics(theName, fsParticles, theDec->motherPart()));
}
else if(theDec->dynType()=="KMatrix"){
std::string pathToConfigFile=theDec->pathToConfigParser();
std::string projectionParticleNames = theDec->projectionParticleNames();
......
//************************************************************************//
// //
// Copyright 2025 Bertram Kopf (bertram@ep1.rub.de) //
// - Ruhr-Universität Bochum //
// //
// This file is part of Pawian. //
// //
// Pawian is free software: you can redistribute it and/or modify //
// it under the terms of the GNU General Public License as published by //
// the Free Software Foundation, either version 3 of the License, or //
// (at your option) any later version. //
// //
// Pawian is distributed in the hope that it will be useful, //
// but WITHOUT ANY WARRANTY; without even the implied warranty of //
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the //
// GNU General Public License for more details. //
// //
// You should have received a copy of the GNU General Public License //
// along with Pawian. If not, see <http://www.gnu.org/licenses/>. //
// //
//************************************************************************//
// SigmaParamDynamics class definition file. -*- C++ -*-
// Copyright 2025 Bertram Kopf
#include <getopt.h>
#include <fstream>
#include <string>
#include "PwaUtils/SigmaParamDynamics.hh"
#include "ErrLogger/ErrLogger.hh"
#include "Particle/Particle.hh"
#include "PwaDynamics/SigmaParameterization.hh"
SigmaParamDynamics::SigmaParamDynamics(std::string& name, std::vector<Particle*>& fsParticles, Particle* mother) :
AbsDynamics(name, fsParticles, mother)
,_currentStrength(1.)
,_sigmaParamDyn(std::shared_ptr<SigmaParameterization>(new SigmaParameterization))
{
_isLdependent=false;
}
SigmaParamDynamics::~SigmaParamDynamics()
{
}
complex<double> SigmaParamDynamics::eval(EvtData* theData, AbsXdecAmp* grandmaAmp, Spin OrbMom){
if ( !_recalculate){
return _cachedMap[theData->evtNo];
}
complex<double> result=_sigmaParamDyn->calc(theData->DoubleMassId.at(_dynId), _currentStrength);
if ( _cacheAmps){
_cachedMap[theData->evtNo]=result;
}
return result;
}
void SigmaParamDynamics::fillDefaultParams(std::shared_ptr<AbsPawianParameters> fitPar){
//fill mass
std::string strengthName=_massKey+"Strength";
fitPar->Add(strengthName, _currentStrength, 0.1);
// fitPar->SetLimits(massName, minMass, maxMass);
}
void SigmaParamDynamics::updateFitParams(std::shared_ptr<AbsPawianParameters> fitPar){
std::string strengthName=_massKey+"Strength";
_currentStrength=fitPar->Value(strengthName);
}
void SigmaParamDynamics::fillParamNameList(){
_paramNameList.clear();
//fill
std::string strengthName=_massKey+"Strength";
_paramNameList.push_back(strengthName);
}
//************************************************************************//
// //
// Copyright 2025 Bertram Kopf (bertram@ep1.rub.de) //
// - Ruhr-Universität Bochum //
// //
// This file is part of Pawian. //
// //
// Pawian is free software: you can redistribute it and/or modify //
// it under the terms of the GNU General Public License as published by //
// the Free Software Foundation, either version 3 of the License, or //
// (at your option) any later version. //
// //
// Pawian is distributed in the hope that it will be useful, //
// but WITHOUT ANY WARRANTY; without even the implied warranty of //
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the //
// GNU General Public License for more details. //
// //
// You should have received a copy of the GNU General Public License //
// along with Pawian. If not, see <http://www.gnu.org/licenses/>. //
// //
//************************************************************************//
// SigmaParamDynamics class definition file. -*- C++ -*-
// Copyright 2025 Bertram Kopf
#pragma once
#include <iostream>
#include <vector>
#include <complex>
#include <map>
#include <string>
#include <memory>
#include "PwaUtils/AbsDynamics.hh"
#include "FitParams/AbsPawianParameters.hh"
class SigmaParameterization;
class SigmaParamDynamics : public AbsDynamics{
public:
SigmaParamDynamics(std::string& name, std::vector<Particle*>& fsParticles, Particle* mother);
virtual ~SigmaParamDynamics();
virtual std::string type() {return "SigmaParamDynamics";}
virtual complex<double> eval(EvtData* theData, AbsXdecAmp* grandmaAmp, Spin OrbMom=0);
virtual void fillDefaultParams(std::shared_ptr<AbsPawianParameters> fitPar);
virtual void updateFitParams(std::shared_ptr<AbsPawianParameters> fitPar);
virtual void fillParamNameList();
protected:
// std::string _massKey;
double _currentStrength;
std::shared_ptr<SigmaParameterization> _sigmaParamDyn;
private:
};
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment