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

added LASS parametrization for KpiS-wave below 1.6 GeV/c2

parent f67a63e4
No related branches found
No related tags found
No related merge requests found
......@@ -295,6 +295,9 @@ void DecAngularDist::fillTensorAmps(double cosT){
Tensor<complex<double> > epsilonDaughter1Project = _polDaughter1(lamDaughter1);
for (Spin lamDaughter2=-_daughter2JPC->J; lamDaughter2<=_daughter2JPC->J; ++lamDaughter2){
// if( fabs(lamMother-lamDaughter1-lamDaughter2)> (*itLS).L) continue;
// if( fabs(lamDaughter1+lamDaughter2)> (*itLS).S) continue;
// Spin lamDaughterDiff=lamDaughter1-lamDaughter2;
// if( fabs(lamDaughterDiff)>_motherJPC->J || fabs(lamDaughterDiff)>(*itLS).S) continue;
......
......@@ -37,6 +37,8 @@
#include "PwaDynamics/KPole.hh"
//#include "PwaDynamics/KPoleNonRel.hh"
#include "PwaDynamics/KMatrixRel.hh"
#include "PwaDynamics/LASS.hh"
#include "TFile.h"
#include "TH1F.h"
......@@ -52,6 +54,17 @@ KPiSWaveTMatrix::KPiSWaveTMatrix() :
std::string rootFileName="./KPiSWaveTMatrix.root";
_theTFile=new TFile(rootFileName.c_str(),"recreate");
double mLASS=1.435;
double gammaLASS=0.279;
// double aLASS=1.9/1000.;
// double rLASS=1.76/1000.;
double aLASS=1.9;
double rLASS=1.76;
double BLASS=1.0;
double phiB=0.;
double RLASS=1.0;
double phiR=0.;
int size=2000;
double massMin=0.135+0.5;
double massMax=1.8;
......@@ -63,6 +76,10 @@ KPiSWaveTMatrix::KPiSWaveTMatrix() :
_KPiAmpImagH1= new TH1F("_KPiAmpImagH1","K #pi amp Im",size+1, massMin, massMax);
_KPiAmpImagH1->SetYTitle("K #pi amp Im");
_KPiAmpRealLASSH1= new TH1F("_KPiAmpRealLASSH1","K #pi amp Re LASS",size+1, massMin, massMax);
_KPiAmpRealLASSH1->SetYTitle("K #pi amp Re");
_KPiAmpImagLASSH1= new TH1F("_KPiAmpImagLASSH1","K #pi amp Im LASS",size+1, massMin, massMax);
_KPiAmpImagLASSH1->SetYTitle("K #pi amp Im");
std::shared_ptr<KMatrixBase> theKMatrixIso12(new KMatrixKPiSFocus(1));
std::shared_ptr<KMatrixBase> theKMatrixIso32(new KMatrixKPiSFocus(3));
......@@ -75,7 +92,11 @@ KPiSWaveTMatrix::KPiSWaveTMatrix() :
for (double mass=massMin; mass<massMax; mass+=stepSize){
Vector4<double> mass4Vec(mass, 0.,0.,0.);
complex<double> currentLASS=LASS::K0star_1430(mass, mLASS, gammaLASS, aLASS, rLASS, BLASS, phiB, RLASS, phiR);
_KPiAmpRealLASSH1->Fill(mass, currentLASS.real());
_KPiAmpImagLASSH1->Fill(mass, currentLASS.imag());
theTMatrix12->evalMatrix(mass);
theTMatrix32->evalMatrix(mass);
......
......@@ -66,6 +66,8 @@ private:
TFile* _theTFile;
TH1F* _KPiAmpRealH1;
TH1F* _KPiAmpImagH1;
TH1F* _KPiAmpRealLASSH1;
TH1F* _KPiAmpImagLASSH1;
};
//************************************************************************//
// //
// Copyright 2014 Bertram Kopf (bertram@ep1.rub.de) //
// Julian Pychy (julian@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/>. //
// //
//************************************************************************//
#include "qft++/relativistic-quantum-mechanics/Utils.hh"
#include "PwaDynamics/LASS.hh"
complex<double> LASS::K0star_1430(double currentMass, double m0, double gammaM, double aLASS, double rLASS, double BLASS, double phiB, double RLASS, double phiR){
complex<double> result(0.,0.);
if(currentMass>1.6) return result;
double KMass=0.49367;
double piMass=0.13498;
double breakupQ=breakupMomQ(currentMass, KMass, piMass).real();
//calculate the background phase motion
double cot_deltaB = 1.0/(aLASS*breakupQ) + 0.5*rLASS*breakupQ;
double deltaB = atan( 1.0/cot_deltaB);
double totalB = deltaB + phiB ;
//calculate the resonant phase motion
double deltaR = atan((m0*gammaM/(m0*m0 - currentMass*currentMass)));
double totalR = deltaR + phiR ;
//sum them up
complex<double> bkgB = complex<double>(BLASS*sin(totalB),0)*complex<double>(cos(totalB),sin(totalB));
complex<double> resT = complex<double>(rLASS*sin(deltaR),0)*complex<double>(cos(totalR),sin(totalR))*complex<double>(cos(2*totalB),sin(2*totalB));
result = bkgB + resT;
return result;
}
//************************************************************************//
// //
// Copyright 2013 Bertram Kopf (bertram@ep1.rub.de) //
// Julian Pychy (julian@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/>. //
// //
//************************************************************************//
// Copyright 2014 Bertram Kopf
#pragma once
//_____________________________________________________________________________
// @file LASS.hh
//_____________________________________________________________________________
#include <iostream>
#include <complex>
#include <utility>
#include <memory>
#include "qft++/relativistic-quantum-mechanics/Utils.hh"
//_____________________________________________________________________________
//_____________________________________________________________________________
namespace LASS {
complex<double> K0star_1430(double currentMass, double m0, double gammaM, double aLASS, double rLASS, double BLASS, double phiB, double RLASS, double phiR);
}; // namespace LASS
......@@ -40,6 +40,8 @@
#include "PwaUtils/PiPiSWaveASDynamics.hh"
#include "PwaUtils/KMatrixDynamics.hh"
#include "PwaUtils/VoigtDynamics.hh"
#include "PwaUtils/K0star1430LassDynamics.hh"
#include "PwaUtils/GlobalEnv.hh"
#include "PwaUtils/WoDynamics.hh"
......@@ -103,6 +105,8 @@ std::shared_ptr<AbsDynamics> DynRegistry::getDynamics(std::shared_ptr<AbsDecay>
result= std::shared_ptr<AbsDynamics>(new PiPiSWaveASDynamics(theName, fsParticles, theDec->motherPart(), GlobalEnv::instance()->particleTable()));
else if(theDec->dynType()=="Voigt")
result= std::shared_ptr<AbsDynamics>(new VoigtDynamics(theName, fsParticles, theDec->motherPart()));
else if(theDec->dynType()=="K0star1430Lass")
result= std::shared_ptr<AbsDynamics>(new K0star1430LassDynamics(theName, fsParticles, theDec->motherPart()));
else if(theDec->dynType()=="WoDynamics") result= std::shared_ptr<AbsDynamics>(new WoDynamics(theName, fsParticles, theDec->motherPart()));
else{
Alert << "Dyn type:\t" << theDec->dynType() << "\tdoes not exist" << endmsg;
......
//************************************************************************//
// //
// Copyright 2013 Bertram Kopf (bertram@ep1.rub.de) //
// Julian Pychy (julian@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/>. //
// //
//************************************************************************//
// K0star1430LassDynamics class definition file. -*- C++ -*-
// Copyright 2013 Bertram Kopf
#include <getopt.h>
#include <fstream>
#include <string>
#include <mutex>
#include <boost/algorithm/string.hpp>
#include "PwaUtils/K0star1430LassDynamics.hh"
#include "ErrLogger/ErrLogger.hh"
#include "Particle/Particle.hh"
#include "PwaDynamics/LASS.hh"
#include "Utils/FunctionUtils.hh"
K0star1430LassDynamics::K0star1430LassDynamics(std::string& name, std::vector<Particle*>& fsParticles, Particle* mother) :
AbsDynamics(name, fsParticles, mother)
,_aLASSKey(_massKey+"aLASS_PosNeg")
,_rLASSKey(_massKey+"rLASS_PosNeg")
,_BLASSKey(_massKey+"BLASS_PosNeg")
,_phiBKey(_massKey+"phiB_PosNeg")
,_RLASSKey(_massKey+"RLASS_PosNeg")
,_phiRKey(_massKey+"phiR_PosNeg")
{
}
K0star1430LassDynamics::~K0star1430LassDynamics()
{
}
complex<double> K0star1430LassDynamics::eval(EvtData* theData, AbsXdecAmp* grandmaAmp, Spin OrbMom){
int evtNo=theData->evtNo;
if ( _cacheAmps && !_recalculate){
return _cachedMap[evtNo];
}
complex<double> result=LASS::K0star_1430(theData->DoubleString.at(_dynKey), _currentMass, _currentWidth, _currentaLASS, _currentrLASS, _currentBLASS, _currentphiB, _currentRLASS, _currentphiR);
if ( _cacheAmps){
theMutex.lock();
_cachedMap[evtNo]=result;
theMutex.unlock();
}
return result;
}
void K0star1430LassDynamics::getDefaultParams(fitParams& fitVal, fitParams& fitErr){
fitVal.Masses[_massKey]=_mother->mass();
fitErr.Masses[_massKey]=0.03;
fitVal.Widths[_massKey]=_mother->width();
fitErr.Widths[_massKey]=0.2*_mother->width();
fitVal.otherParams[_aLASSKey]=1.07;
fitErr.otherParams[_aLASSKey]=0.05;
fitVal.otherParams[_rLASSKey]=-2.852;
fitErr.otherParams[_rLASSKey]=0.05;
fitVal.otherParams[_BLASSKey]=0.7;
fitErr.otherParams[_BLASSKey]=0.05;
fitVal.otherParams[_phiBKey]=0.7;
fitErr.otherParams[_phiBKey]=0.05;
fitVal.otherParams[_RLASSKey]=1.;
fitErr.otherParams[_RLASSKey]=0.05;
fitVal.otherParams[_phiRKey]=-5.356;
fitErr.otherParams[_phiRKey]=0.05;
}
bool K0star1430LassDynamics::checkRecalculation(fitParams& theParamVal){
_recalculate=false;
double mass=theParamVal.Masses[_massKey];
if (!CheckDoubleEquality(mass, _currentMass)){
_recalculate=true;
return _recalculate;
}
double width=theParamVal.Widths[_massKey];
if (!CheckDoubleEquality(width, _currentWidth)){
_recalculate=true;
return _recalculate;
}
double aLASS=theParamVal.otherParams[_aLASSKey];
if (!CheckDoubleEquality(aLASS, _currentaLASS)){
_recalculate=true;
return _recalculate;
}
double rLASS=theParamVal.otherParams[_rLASSKey];
if (!CheckDoubleEquality(rLASS, _currentrLASS)){
_recalculate=true;
return _recalculate;
}
double BLASS=theParamVal.otherParams[_BLASSKey];
if (!CheckDoubleEquality(BLASS, _currentBLASS)){
_recalculate=true;
return _recalculate;
}
double phiB=theParamVal.otherParams[_phiBKey];
if (!CheckDoubleEquality(phiB, _currentphiB)){
_recalculate=true;
return _recalculate;
}
double RLASS=theParamVal.otherParams[_RLASSKey];
if (!CheckDoubleEquality(RLASS, _currentRLASS)){
_recalculate=true;
return _recalculate;
}
double phiR=theParamVal.otherParams[_phiRKey];
if (!CheckDoubleEquality(phiR, _currentphiR)){
_recalculate=true;
return _recalculate;
}
return _recalculate;
}
void K0star1430LassDynamics::updateFitParams(fitParams& theParamVal){
_currentMass=theParamVal.Masses.at(_massKey);
_currentWidth=theParamVal.Widths.at(_massKey);
_currentaLASS=theParamVal.otherParams.at(_aLASSKey);
_currentrLASS=theParamVal.otherParams.at(_rLASSKey);
_currentBLASS=theParamVal.otherParams.at(_BLASSKey);
_currentphiB=theParamVal.otherParams.at(_phiBKey);
_currentRLASS=theParamVal.otherParams.at(_RLASSKey);
_currentphiR=theParamVal.otherParams.at(_phiRKey);
}
void K0star1430LassDynamics::setMassKey(std::string& theMassKey){
boost::replace_all(_aLASSKey,_massKey, theMassKey);
boost::replace_all(_rLASSKey,_massKey, theMassKey);
boost::replace_all(_BLASSKey,_massKey, theMassKey);
boost::replace_all(_phiBKey,_massKey, theMassKey);
boost::replace_all(_RLASSKey,_massKey, theMassKey);
boost::replace_all(_phiRKey,_massKey, theMassKey);
AbsDynamics::setMassKey(theMassKey);
}
//************************************************************************//
// //
// Copyright 2013 Bertram Kopf (bertram@ep1.rub.de) //
// Julian Pychy (julian@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/>. //
// //
//************************************************************************//
// K0star1430LassDynamics class definition file. -*- C++ -*-
// Copyright 2014 Bertram Kopf
#pragma once
#include <iostream>
#include <vector>
#include <complex>
#include <map>
#include <string>
#include <memory>
#include "PwaUtils/AbsDynamics.hh"
class K0star1430Lass;
class K0star1430LassDynamics : public AbsDynamics{
public:
K0star1430LassDynamics(std::string& name, std::vector<Particle*>& fsParticles, Particle* mother);
virtual ~K0star1430LassDynamics();
virtual complex<double> eval(EvtData* theData, AbsXdecAmp* grandmaAmp, Spin OrbMom=0);
virtual void getDefaultParams(fitParams& fitVal, fitParams& fitErr);
virtual bool checkRecalculation(fitParams& theParamVal);
virtual void updateFitParams(fitParams& theParamVal);
virtual void setMassKey(std::string& theMassKey);
protected:
// std::string _key;
std::string _aLASSKey;
std::string _rLASSKey;
std::string _BLASSKey;
std::string _phiBKey;
std::string _RLASSKey;
std::string _phiRKey;
double _currentMass;
double _currentWidth;
double _currentaLASS;
double _currentrLASS;
double _currentBLASS;
double _currentphiB;
double _currentRLASS;
double _currentphiR;
std::map<int, complex<double> > _cachedMap;
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