else if (TString(name) == "CuHN") {fN = 63; fR = 4.28; fA = 0.5; fW = 0; fF = 1; fZ=29;} // from arXiv:0904.4080v1
else if (TString(name) == "Xe") {fN = 129; fR = 5.36; fA = 0.59; fW = 0; fF = 1; fZ=54;} // adapted from arXiv:1703.04278
else if (TString(name) == "Xes") {fN = 129; fR = 5.42; fA = 0.57; fW = 0; fF = 1; fZ=54;} // scale from Sb (Antimony, A=122, r=5.32) by 1.019 = (129/122)**0.333
else if (TString(name) == "Xe2") {fN = 129; fR = 5.36; fA = 0.59; fW = 0; fF = 8; fZ=54; fBeta2=0.161; fBeta4=-0.003;} // adapted from arXiv:1703.04278 and Z. Physik (1974) 270: 113
else if (TString(name) == "Xe2a") {fN = 129; fR = 5.36; fA = 0.59; fW = 0; fF = 8; fZ=54; fBeta2=0.18; fBeta4=0;} // ALICE parameters (see public note from 2018 at https://cds.cern.ch/collection/ALICE%20Public%20Notes?ln=en)
else if (TString(name) == "Pbrw") {fN = 208; fR = 6.62; fA = 0.546; fW = 0; fF = 12; fZ=82; r0=1.00863; r1=-0.00044808; r2=-0.000205872;} //only Pb 207 was tested but should be the same for 208
cout << "Could not find nucleus " << name << endl;
return;
}
switch (fF) {
case 0: // Proton exp
fFunc1 = new TF1("prot","x*x*exp(-x/[0])",0,5);
fFunc1->SetParameter(0,fR);
break;
case 1: // 3pF
fFunc1 = new TF1(name,"x*x*(1+[2]*(x/[0])**2)/(1+exp((x-[0])/[1]))",0,fMaxR);
fFunc1->SetParameters(fR,fA,fW);
break;
case 2: // 3pG
fFunc1 = new TF1("3pg","x*x*(1+[2]*(x/[0])**2)/(1+exp((x**2-[0]**2)/[1]**2))",0,fMaxR);
fFunc1->SetParameters(fR,fA,fW);
break;
case 3: // Hulthen (see nucl-ex/0603010)
fFunc1 = new TF1("f3","x*x*([0]*[1]*([0]+[1]))/(2*pi*(pow([0]-[1],2)))*pow((exp(-[0]*x)-exp(-[1]*x))/x,2)",0,fMaxR);
fFunc1->SetParameters(1/4.38,1/.85);
break;
case 4: // Hulthen HIJING
fFunc1 = new TF1("f4","x*x*([0]*[1]*([0]+[1]))/(2*pi*(pow([0]-[1],2)))*pow((exp(-[0]*x)-exp(-[1]*x))/x,2)",0,fMaxR);
fFunc1->SetParameters(2/4.38,2/.85);
break;
case 5: // Ellipsoid (Uranium)
fFunc1 = new TF1(name,"x*x*(1+[2]*(x/[0])**2)/(1+exp((x-[0])/[1]))",0,fMaxR);
fFunc1->SetParameters(fR,fA,0); // same as 3pF but setting W to zero
break;
case 6: // He3/H3
fFunc1 = 0; // read in file instead
break;
case 7: // Deformed nuclei, box method
#ifndef HAVE_MATHMORE
cerr << "Need libMathMore.so for deformed nuclei" << endl;
gSystem->Exit(123);
#endif
fFunc1 = 0; // no func: only need beta parameters and use uniform box distribution
break;
case 8: // Deformed nuclei, TF2 method
fFunc3 = new TF2("f77","x*x*TMath::Sin(y)/(1+exp((x-[0]*(1+[2]*0.315*(3*pow(cos(y),2)-1.0)+[3]*0.105*(35*pow(cos(y),4)-30*pow(cos(y),2)+3)))/[1]))",0,fMaxR,0.0,TMath::Pi());
fFunc3->SetNpx(120);
fFunc3->SetNpy(120);
fFunc3->SetParameters(fR,fA,fBeta2,fBeta4);
break;
case 9: // Proton gaus
fFunc1 = new TF1("prot","x*x*exp(-x*x/[0]/[0]/2)",0,5);
fFunc1->SetParameter(0,fR);
break;
case 10: // Proton dgaus
fFunc1 = new TF1("prot","x*x*((1-[0])/[1]^3*exp(-x*x/[1]/[1])+[0]/(0.4*[1])^3*exp(-x*x/(0.4*[1])^2))",0,5);
fFunc1->SetParameter(0,0.5);
fFunc1->SetParameter(1,fR);
break;
case 11: // 3pF for proton and neutrons
fFunc1 = new TF1(name,"x*x*(1+[2]*(x/[0])**2)/(1+exp((x-[0])/[1]))",0,fMaxR);
fFunc1->SetParameters(fR,fA,fW);
fFunc2 = new TF1(name,"x*x*(1+[2]*(x/[0])**2)/(1+exp((x-[0])/[1]))",0,fMaxR);
fFunc2->SetParameters(fR2,fA2,fW2);
break;
case 12: // reweighted
fFunc1 = new TF1(name,"x*x*(1+[2]*(x/[0])**2)/(1+exp((x-[0])/[1]))/([3]+[4]*x+[5]*x^2)",0,fMaxR);
fFunc1->SetParameters(fR,fA,fW,r0,r1,r2);
fRecenter=1;
fSmax=0.1;
break;
case 13: // Pb for proton and neutrons reweighted
fFunc1 = new TF1(Form("%s_prot",name),"x*x*(1+[2]*(x/[0])**2)/(1+exp((x-[0])/[1]))/([3]+[4]*x+[5]*x^2)",0,fMaxR);
case 14: // Deformed nuclei, TF2 method, reweighted
fFunc3 = new TF2("f77","x*x*TMath::Sin(y)/(1+exp((x-[0]*(1+[2]*0.315*(3*pow(cos(y),2)-1.0)+[3]*0.105*(35*pow(cos(y),4)-30*pow(cos(y),2)+3)))/[1]))/([4]+[5]*x+[6]*x^2)",0,fMaxR,0.0,TMath::Pi());