Commit a48fc554 authored by Frank Siegert's avatar Frank Siegert

fixes and temporary hacks to make (partonic) hadron decays work again on trunk.

parent 1595cf59
......@@ -58,6 +58,7 @@ Particle_Info::Particle_Info
m_formfactor(0), m_on(on), m_massive(1), m_hadron(1), m_isgroup(0),
m_idname(idname), m_antiname(antiname)
{
m_antiname=m_idname+"b";
m_content.push_back(new Flavour(*this));
}
......@@ -70,6 +71,7 @@ Particle_Info::Particle_Info
m_formfactor(formfactor), m_on(1), m_massive(1), m_hadron(1), m_isgroup(0),
m_idname(idname), m_antiname(antiname)
{
m_antiname=m_idname+"b";
m_content.push_back(new Flavour(*this));
}
......
......@@ -387,7 +387,7 @@ void Beam_Spectra_Handler::InitializeFlav(kf_code flav)
}
else if (flav==kf_n) {
s_kftable[flav]=
new Particle_Info(kf_n,0.939566,7.424e-28,0,0,1,1,"n","n");
new Particle_Info(kf_n,0.939566,7.424e-28,0,0,1,1,"n","nb");
initialize_diquarks=true;
}
else if (flav==kf_e) {
......@@ -425,31 +425,31 @@ void Beam_Spectra_Handler::InitializeFlav(kf_code flav)
if (initialize_diquarks) {
s_kftable[1103]=new Particle_Info(1103,0.77133,0,-2,-3,2,0,0,1,0,"dd_1","dd_1");
s_kftable[2101]=new Particle_Info(2101,0.57933,0,1,-3,0,0,0,1,0,"ud_0","ud_0");
s_kftable[2103]=new Particle_Info(2103,0.77133,0,1,-3,2,0,0,1,0,"ud_1","ud_1");
s_kftable[2203]=new Particle_Info(2203,0.77133,0,4,-3,2,0,0,1,0,"uu_1","uu_1");
s_kftable[3101]=new Particle_Info(3101,0.80473,0,-2,-3,0,0,0,1,0,"sd_0","sd_0");
s_kftable[3103]=new Particle_Info(3103,0.92953,0,-2,-3,2,0,0,1,0,"sd_1","sd_1");
s_kftable[3201]=new Particle_Info(3201,0.80473,0,1,-3,0,0,0,1,0,"su_0","su_0");
s_kftable[3203]=new Particle_Info(3203,0.92953,0,1,-3,2,0,0,1,0,"su_1","su_1");
s_kftable[3303]=new Particle_Info(3303,1.09361,0,-2,-3,2,0,0,1,0,"ss_1","ss_1");
s_kftable[4101]=new Particle_Info(4101,1.96908,0,1,-3,0,0,0,1,0,"cd_0","cd_0");
s_kftable[4103]=new Particle_Info(4103,2.00808,0,1,-3,2,0,0,1,0,"cd_1","cd_1");
s_kftable[4201]=new Particle_Info(4201,1.96908,0,4,-3,0,0,0,1,0,"cu_0","cu_0");
s_kftable[4203]=new Particle_Info(4203,2.00808,0,4,-3,2,0,0,1,0,"cu_1","cu_1");
s_kftable[4301]=new Particle_Info(4301,2.15432,0,1,-3,0,0,0,1,0,"cs_0","cs_0");
s_kftable[4303]=new Particle_Info(4303,2.17967,0,1,-3,2,0,0,1,0,"cs_1","cs_1");
s_kftable[4403]=new Particle_Info(4403,3.27531,0,4,-3,2,0,0,1,0,"cc_1","cc_1");
s_kftable[5101]=new Particle_Info(5101,5.38897,0,-2,-3,0,0,0,1,0,"bd_0","bd_0");
s_kftable[5103]=new Particle_Info(5103,5.40145,0,-2,-3,2,0,0,1,0,"bd_1","bd_1");
s_kftable[5201]=new Particle_Info(5201,5.38897,0,1,-3,0,0,0,1,0,"bu_0","bu_0");
s_kftable[5203]=new Particle_Info(5203,5.40145,0,1,-3,2,0,0,1,0,"bu_1","bu_1");
s_kftable[5301]=new Particle_Info(5301,5.56725,0,-2,-3,0,0,0,1,0,"bs_0","bs_0");
s_kftable[5303]=new Particle_Info(5303,5.57536,0,-2,-3,2,0,0,1,0,"bs_1","bs_1");
s_kftable[5401]=new Particle_Info(5401,6.67143,0,1,-3,0,0,0,1,0,"bc_0","bc_0");
s_kftable[5403]=new Particle_Info(5403,6.67397,0,1,-3,2,0,0,1,0,"bc_1","bc_1");
s_kftable[5503]=new Particle_Info(5503,10.07354,0,-2,-3,2,0,0,1,0,"bb_1","bb_1");
s_kftable[1103]=new Particle_Info(1103,0.77133,0,-2,-3,2,0,0,1,0,"dd_1","dd_1b");
s_kftable[2101]=new Particle_Info(2101,0.57933,0,1,-3,0,0,0,1,0,"ud_0","ud_0b");
s_kftable[2103]=new Particle_Info(2103,0.77133,0,1,-3,2,0,0,1,0,"ud_1","ud_1b");
s_kftable[2203]=new Particle_Info(2203,0.77133,0,4,-3,2,0,0,1,0,"uu_1","uu_1b");
s_kftable[3101]=new Particle_Info(3101,0.80473,0,-2,-3,0,0,0,1,0,"sd_0","sd_0b");
s_kftable[3103]=new Particle_Info(3103,0.92953,0,-2,-3,2,0,0,1,0,"sd_1","sd_1b");
s_kftable[3201]=new Particle_Info(3201,0.80473,0,1,-3,0,0,0,1,0,"su_0","su_0b");
s_kftable[3203]=new Particle_Info(3203,0.92953,0,1,-3,2,0,0,1,0,"su_1","su_1b");
s_kftable[3303]=new Particle_Info(3303,1.09361,0,-2,-3,2,0,0,1,0,"ss_1","ss_1b");
s_kftable[4101]=new Particle_Info(4101,1.96908,0,1,-3,0,0,0,1,0,"cd_0","cd_0b");
s_kftable[4103]=new Particle_Info(4103,2.00808,0,1,-3,2,0,0,1,0,"cd_1","cd_1b");
s_kftable[4201]=new Particle_Info(4201,1.96908,0,4,-3,0,0,0,1,0,"cu_0","cu_0b");
s_kftable[4203]=new Particle_Info(4203,2.00808,0,4,-3,2,0,0,1,0,"cu_1","cu_1b");
s_kftable[4301]=new Particle_Info(4301,2.15432,0,1,-3,0,0,0,1,0,"cs_0","cs_0b");
s_kftable[4303]=new Particle_Info(4303,2.17967,0,1,-3,2,0,0,1,0,"cs_1","cs_1b");
s_kftable[4403]=new Particle_Info(4403,3.27531,0,4,-3,2,0,0,1,0,"cc_1","cc_1b");
s_kftable[5101]=new Particle_Info(5101,5.38897,0,-2,-3,0,0,0,1,0,"bd_0","bd_0b");
s_kftable[5103]=new Particle_Info(5103,5.40145,0,-2,-3,2,0,0,1,0,"bd_1","bd_1b");
s_kftable[5201]=new Particle_Info(5201,5.38897,0,1,-3,0,0,0,1,0,"bu_0","bu_0b");
s_kftable[5203]=new Particle_Info(5203,5.40145,0,1,-3,2,0,0,1,0,"bu_1","bu_1b");
s_kftable[5301]=new Particle_Info(5301,5.56725,0,-2,-3,0,0,0,1,0,"bs_0","bs_0b");
s_kftable[5303]=new Particle_Info(5303,5.57536,0,-2,-3,2,0,0,1,0,"bs_1","bs_1b");
s_kftable[5401]=new Particle_Info(5401,6.67143,0,1,-3,0,0,0,1,0,"bc_0","bc_0b");
s_kftable[5403]=new Particle_Info(5403,6.67397,0,1,-3,2,0,0,1,0,"bc_1","bc_1b");
s_kftable[5503]=new Particle_Info(5503,10.07354,0,-2,-3,2,0,0,1,0,"bb_1","bb_1b");
}
}
}
......
......@@ -38,35 +38,29 @@ void Hadron_Decay_Table::Read(std::string path, std::string file)
abort();
}
int nchannels(0), rewrite(0);
vector<int> helpkfc;
double BR, dBR, totBR(0.);
string origin;
bool haspartonics(false);
std::vector<int> specs;
std::vector<double> specweights;
for (size_t i=0;i<helpsvv.size();i++) {
if ( helpsvv[i][0] == string("NO_ANTI") )
continue;
if (helpsvv[i].size()!=3 && helpsvv[i][0]==string("SPECTATORS")) {
haspartonics = true;
Tools::ExtractSpecs(helpsvv[i][1],specs,specweights);
msg_Tracking()<<METHOD<<": found spectators for "<<m_flin<<": ";
for (size_t j=0;j<specs.size();j++) msg_Tracking()<<specs[j]<<" ";
msg_Tracking()<<"\n";
DEBUG_INFO("found spectators for "<<m_flin);
for (size_t j=0;j<specs.size();j++) DEBUG_VAR(specs[j]);
}
}
int nchannels(0), rewrite(0);
vector<int> helpkfc;
double BR, dBR, totBR(0.);
string origin;
Flavour flav;
Hadron_Decay_Channel* hdc;
for (size_t i=0;i<helpsvv.size();i++) {
if ( helpsvv[i][0] == string("NO_ANTI") )
continue;
if (helpsvv[i].size()>1 && Tools::ExtractFlavours(helpkfc,helpsvv[i][0])) {
else if (helpsvv[i].size()>1 && Tools::ExtractFlavours(helpkfc,helpsvv[i][0])) {
Tools::ExtractBRInfo(helpsvv[i][1], BR, dBR, origin);
hdc = new Hadron_Decay_Channel(Flav(),p_ms,path);
Hadron_Decay_Channel* hdc = new Hadron_Decay_Channel(Flav(),p_ms,path);
int charge = Flav().IntCharge();
double mass = Flav().HadMass();
for (size_t j=0;j<helpkfc.size();++j) {
flav = Flavour(abs(helpkfc[j]));
Flavour flav = Flavour(abs(helpkfc[j]));
if (helpkfc[j]<0) flav = flav.Bar();
hdc->AddDecayProduct(flav);
charge-=flav.IntCharge();
......@@ -76,7 +70,7 @@ void Hadron_Decay_Table::Read(std::string path, std::string file)
msg_Tracking()<<"Found too low mass.";
BR = 0.; dBR = 0.; continue;
}
if (haspartonics) totBR += BR;
totBR += BR;
hdc->SetWidth(BR*Flav().Width());
hdc->SetDeltaWidth(dBR*Flav().Width());
hdc->SetOrigin(origin);
......@@ -94,9 +88,8 @@ void Hadron_Decay_Table::Read(std::string path, std::string file)
nchannels++;
}
}
if (Flav().IsC_Hadron())
msg_Tracking()<<METHOD<<" for "<<om::green<<Flav()<<om::reset<<": "
<<"total BR = "<<om::green<<totBR<<om::reset<<".\n";
DEBUG_VAR(totBR);
if (haspartonics) {
PHASIC::Decay_Table * dectable(NULL);
if (!Flav().IsB_Hadron() && !Flav().IsC_Hadron()) {
......@@ -104,24 +97,14 @@ void Hadron_Decay_Table::Read(std::string path, std::string file)
<<" No suitable partonic decay table found for "
<<Flav()<<".\n"
<<" Will continue and hope for the best.\n";
if(rewrite) {
PRINT_INFO("TODO: migrate to Decaydata.db");
/*
Move(path+file, path+"."+file+".old");
ofstream ostr( (path + file).c_str() );
Write(ostr);
ostr.close();
*/
}
ScaleToWidth();
return;
}
double totspec(0.);
Flavour spec;
for (size_t k=0;k<specs.size();k++) totspec+=specweights[k];
for (size_t k=0;k<specs.size();k++) {
bool isAnti(false);
spec = Flavour(abs(specs[k]));
Flavour spec = Flavour(abs(specs[k]));
if (specs[k]<0) spec = spec.Bar();
if ((spec.IsQuark() && !spec.IsAnti()) ||
(spec.IsDiQuark() && spec.IsAnti())) isAnti=true;
......@@ -160,38 +143,28 @@ void Hadron_Decay_Table::Read(std::string path, std::string file)
msg_Tracking()<<".\n";
continue;
}
hdc = new Hadron_Decay_Channel(Flav(),p_ms,path);
Hadron_Decay_Channel* hdc = new Hadron_Decay_Channel(Flav(),p_ms,path);
hdc->AddDecayProduct(spec,false);
std::string filename=m_flin.IDName();
if (filename.find("_{s}")!=string::npos)
filename = StringReplace(filename, "_{s}", "s");
if (filename.find("-_{b}")!=string::npos)
filename = StringReplace(filename, "-_{b}", "b-");
if (filename.find("_{b}")!=string::npos)
filename = StringReplace(filename, "_{b}", "b");
if (filename.find("++_{c}")!=string::npos)
filename = StringReplace(filename, "++_{c}", "c++");
if (filename.find("+_{c}")!=string::npos)
filename = StringReplace(filename, "+_{c}", "c+");
if (filename.find("_{c}")!=string::npos)
filename = StringReplace(filename, "_{c}", "c");
if (filename.find("_fict")!=string::npos)
filename = StringReplace(filename, "_fict", "");
filename += "_"+spec.IDName();
std::string filename=m_flin.ShellName();
filename += "_"+spec.ShellName();
msg_Tracking()<<" Add partonic decay: "<<Flav()<<" --> ";
for (size_t j=0;j<(*dectable)[i]->NOut();j++) {
flav = (*dectable)[i]->GetDecayProduct(j);
Flavour flav = (*dectable)[i]->GetDecayProduct(j);
if (isAnti) flav=flav.Bar();
msg_Tracking()<<flav<<" ";
hdc->AddDecayProduct(flav,false);
charge -= flav.IntCharge();
mass -= flav.HadMass();
filename += flav.IDName();
filename += flav.ShellName();
}
hdc->SetWidth(BR*partWidth);
hdc->SetDeltaWidth(0.);
hdc->SetOrigin("");
filename += ".dat";
// temporary hack to avoid filename changes in Decaydata.db
filename = StringReplace(filename, "ve", "nu_e");
filename = StringReplace(filename, "vmu", "nu_mu");
filename = StringReplace(filename, "vtau", "nu_tau");
msg_Tracking()<<" ---> "<<filename<<".\n";
hdc->SetFileName(filename);
AddDecayChannel(hdc);
......@@ -208,9 +181,6 @@ void Hadron_Decay_Table::Read(std::string path, std::string file)
ostr.Close();
}
ScaleToWidth();
if (Flav().IsC_Hadron())
msg_Tracking()<<" --> after rescaling: "
<<om::green<<totBR<<om::reset<<".\n";
}
......
......@@ -67,7 +67,6 @@ double Tools::OffShellMassWidth( double s, double Mass2, double Width, double ms
bool Tools::ExtractFlavours(std::vector<int> & helpkfc,std::string help)
{
if (help==std::string("SPECTATORS")) return false;
helpkfc.clear();
size_t pos = help.find("{");
bool hit;
......
......@@ -52,7 +52,7 @@ namespace METOOLS
m_spins = std::vector<int>(flavs.size());
m_n=1;
for(size_t i=0;i<flavs.size();i++) {
if(flavs[i].IsVector() && flavs[i].IsMassive()==0) m_spins[i] = 2;
if(flavs[i].IsVector() && flavs[i].HadMass()==0) m_spins[i] = 2;
else m_spins[i] = flavs[i].IntSpin()+1;
m_n*=m_spins[i];
}
......@@ -66,7 +66,7 @@ namespace METOOLS
m_n=1;
for(size_t i=0;i<indices.size();i++) {
if(flavs[indices[i]].IsVector() && flavs[indices[i]].IsMassive()==0)
if(flavs[indices[i]].IsVector() && flavs[indices[i]].HadMass()==0)
m_spins[i] = 2;
else m_spins[i] = flavs[indices[i]].IntSpin()+1;
m_n*=m_spins[i];
......@@ -80,7 +80,7 @@ namespace METOOLS
m_n=1;
for(size_t i=0;i<particles.size();i++) {
if(particles[i]->Flav().IsVector() && particles[i]->Flav().IsMassive()==0) m_spins[i] = 2;
if(particles[i]->Flav().IsVector() && particles[i]->Flav().HadMass()==0) m_spins[i] = 2;
else m_spins[i] = particles[i]->Flav().IntSpin()+1;
m_n*=m_spins[i];
}
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment