No OneTemporary
Actions

Size

59 KB

Subscribers

None

View Options

	diff --git a/Models/Feynrules/python/ufo2peg/check_lorentz.py b/Models/Feynrules/python/ufo2peg/check_lorentz.py
	--- a/Models/Feynrules/python/ufo2peg/check_lorentz.py
	+++ b/Models/Feynrules/python/ufo2peg/check_lorentz.py
	@@ -1,621 +1,695 @@
	import itertools,cmath,re
	from .helpers import SkipThisVertex
	+from .converter import py2cpp

	# ordering for EW VVV vertices
	def VVVordering(vertex) :
	pattern = "if((p1->id()==%s&&p2->id()==%s&&p3->id()==%s)"+\
	"\|\|(p1->id()==%s&&p2->id()==%s&&p3->id()==%s)\|\|"+\
	"(p1->id()==%s&&p2->id()==%s&&p3->id()==%s)) {norm(-norm());}"
	ordering = pattern%(vertex.particles[1].pdg_code,
	vertex.particles[0].pdg_code,
	vertex.particles[2].pdg_code,
	vertex.particles[0].pdg_code,
	vertex.particles[2].pdg_code,
	vertex.particles[1].pdg_code,
	vertex.particles[2].pdg_code,
	vertex.particles[1].pdg_code,
	vertex.particles[0].pdg_code)
	return ordering

	def tensorCouplings(vertex,value,prefactors,L,lorentztag,pos,all_couplings) :
	# split the structure into its different terms for analysis
	ordering=""
	structure1 = L.structure.split()
	structures =[]
	sign=''
	for struct in structure1 :
	if(struct=='+') :
	continue
	elif(struct=='-') :
	sign='-'
	else :
	structures.append(sign+struct.strip())
	sign=''
	if(lorentztag == 'SST') :
	terms=[['P(1003,2)','P(2003,1)'],
	['P(1003,1)','P(2003,2)'],
	['P(-1,1)','P(-1,2)','Metric(1003,2003)'],
	['Metric(1003,2003)']]
	signs=[1.,1.,-1.,-1.]
	new_couplings=[False]*len(terms)
	elif(lorentztag == 'FFT' ) :
	terms=[['P(2003,1)','Gamma(1003,2,1)'],
	['P(2003,2)','Gamma(1003,2,1)'],
	['P(1003,1)','Gamma(2003,2,1)'],
	['P(1003,2)','Gamma(2003,2,1)'],
	['P(-1,1)','Gamma(-1,2,1)','Metric(1003,2003)'],
	['P(-1,2)','Gamma(-1,2,1)','Metric(1003,2003)'],
	['Metric(1003,2003)']]
	signs=[1.,-1.,1.,-1.,-0.5,0.5,1.]
	new_couplings=[False]3len(terms)
	elif(lorentztag == 'VVT' ) :
	terms=[['P(-1,1)','P(-1,2)','Metric(1,2003)','Metric(2,1003)'], # from C term
	['P(-1,1)','P(-1,2)','Metric(1,1003)','Metric(2,2003)'], # from C term
	['P(-1,1)','P(-1,2)','Metric(1,2)','Metric(1003,2003)'], # from C term
	['P(1,2)','P(2,1)','Metric(1003,2003)'], # from D term (sym)
	['P(1,2)','P(2003,1)','Metric(2,1003)'], # 1st term
	['P(1,2)','P(1003,1)','Metric(2,2003)'], # 1st swap
	['P(2,1)','P(2003,2)','Metric(1,1003)'], # 2nd term
	['P(2,1)','P(1003,2)','Metric(1,2003)'], # 2nd swap
	['P(1003,2)','P(2003,1)','Metric(1,2)'], # 3rd term
	['P(1003,1)','P(2003,2)','Metric(1,2)'], # 3rd swap
	['Metric(1,2003)','Metric(2,1003)'], # from mass term
	['Metric(1,1003)','Metric(2,2003)'], # from mass term
	['Metric(1,2)','Metric(1003,2003)'], # from mass term
	['P(1,1)','P(2,1)','Metric(1003,2003)'], # gauge terms
	['P(1,2)','P(2,2)','Metric(1003,2003)'], # gauge terms
	['P(1,1)','P(2,2)','Metric(1003,2003)'], # gauge terms
	['P(1003,1)','P(1,1)','Metric(2,2003)'], # gauge terms
	['P(1003,2)','P(2,2)','Metric(1,2003)'], # gauge terms
	['P(2003,1)','P(1,1)','Metric(2,1003)'], # gauge terms
	['P(2003,2)','P(2,2)','Metric(1,1003)'], # gauge terms
	]
	signs=[1.,1.,-1.,1.,-1.,-1.,-1.,-1.,1.,1.,1.,1.,-1.,1.,1.,0.25,-1.,-1.,-1.,-1.]
	new_couplings=[False]*len(terms)
	elif(lorentztag == 'FFVT' ) :
	terms = [['Gamma(2004,2,1)','Metric(3,1004)'],
	['Gamma(1004,2,1)','Metric(3,2004)'],
	['Gamma(3,2,1)','Metric(1004,2004)'],
	['Gamma(2004,2,-1)','Metric(3,1004)'],
	['Gamma(1004,2,-1)','Metric(3,2004)'],
	['Gamma(3,2,-1)','Metric(1004,2004)']]
	signs=[1.,1.,-0.5,1.,1.,-0.5]
	new_couplings=[False]3len(terms)
	elif(lorentztag == 'VVVT' ) :
	# the F(mu nu,rho sigma lambda) terms first
	terms = [['P(2004,2)','Metric(1,1004)','Metric(2,3)'],['P(2004,3)','Metric(1,1004)','Metric(2,3)'],
	['P(1004,2)','Metric(1,2004)','Metric(2,3)'],['P(1004,3)','Metric(1,2004)','Metric(2,3)'],
	['P(2004,3)','Metric(1,3)','Metric(2,1004)'],['P(2004,1)','Metric(1,3)','Metric(2,1004)'],
	['P(1004,3)','Metric(1,3)','Metric(2,2004)'],['P(1004,1)','Metric(1,3)','Metric(2,2004)'],
	['P(2004,1)','Metric(1,2)','Metric(3,1004)'],['P(2004,2)','Metric(1,2)','Metric(3,1004)'],
	['P(1004,1)','Metric(1,2)','Metric(3,2004)'],['P(1004,2)','Metric(1,2)','Metric(3,2004)'],
	['P(3,1)','Metric(1,2004)','Metric(2,1004)'],['P(3,2)','Metric(1,2004)','Metric(2,1004)'],
	['P(3,1)','Metric(1,1004)','Metric(2,2004)'],['P(3,2)','Metric(1,1004)','Metric(2,2004)'],
	['P(3,1)','Metric(1,2)','Metric(1004,2004)'],['P(3,2)','Metric(1,2)','Metric(1004,2004)'],
	['P(2,3)','Metric(1,2004)','Metric(3,1004)'],['P(2,1)','Metric(1,2004)','Metric(3,1004)'],
	['P(2,3)','Metric(1,1004)','Metric(3,2004)'],['P(2,1)','Metric(1,1004)','Metric(3,2004)'],
	['P(2,3)','Metric(1,3)','Metric(1004,2004)'],['P(2,1)','Metric(1,3)','Metric(1004,2004)'],
	['P(1,2)','Metric(2,2004)','Metric(3,1004)'],['P(1,3)','Metric(2,2004)','Metric(3,1004)'],
	['P(1,2)','Metric(2,1004)','Metric(3,2004)'],['P(1,3)','Metric(2,1004)','Metric(3,2004)'],
	['P(1,2)','Metric(2,3)','Metric(1004,2004)'],['P(1,3)','Metric(2,3)','Metric(1004,2004)']]
	signs = [1.,-1.,1.,-1.,1.,-1.,1.,-1.,1.,-1.,1.,-1.,
	1.,-1.,1.,-1.,-1.,1.,1.,-1.,1.,-1.,-1.,1.,1.,-1.,1.,-1.,-1.,1.]
	new_couplings=[False]*len(terms)
	l = lambda c: len(pos[c])
	if l(8)!=3 :
	ordering = VVVordering(vertex)
	# unknown
	else :
	raise Exception('Unknown data type "%s".' % lorentztag)
	iterm=0
	for term in terms:
	for perm in itertools.permutations(term):
	label = '*'.join(perm)
	for istruct in range(0,len(structures)) :
	if label in structures[istruct] :
	reminder = structures[istruct].replace(label,'1.',1)
	loc=iterm
	if(reminder.find("ProjM")>=0) :
	reminder=re.sub("\ProjM\(.,.\)","",reminder)
	loc+=len(terms)
	elif(reminder.find("ProjP")>=0) :
	reminder=re.sub("\ProjP\(.,.\)","",reminder)
	loc+=2*len(terms)
	structures[istruct] = "Done"
	value = eval(reminder, {'cmath':cmath} )*signs[iterm]
	if(new_couplings[loc]) :
	new_couplings[loc] += value
	else :
	new_couplings[loc] = value
	iterm+=1
	# check we've handled all the terms
	for val in structures:
	if(val!="Done") :
	raise SkipThisVertex()
	# special for FFVT
	if(lorentztag=="FFVT") :
	t_couplings=new_couplings
	new_couplings=[False]*9
	for i in range(0,9) :
	j = i+3*(i/3)
	k = i+3+3*(i/3)
	if( not t_couplings[j]) :
	new_couplings[i] = t_couplings[k]
	else :
	new_couplings[i] = t_couplings[j]
	# set the couplings
	for icoup in range(0,len(new_couplings)) :
	if(new_couplings[icoup]) :
	new_couplings[icoup] = '(%s) * (%s) *(%s)' % (new_couplings[icoup],prefactors,value)
	if(len(all_couplings)==0) :
	all_couplings=new_couplings
	else :
	for icoup in range(0,len(new_couplings)) :
	if(new_couplings[icoup] and all_couplings[icoup]) :
	all_couplings[icoup] = '(%s) * (%s) *(%s) + (%s) ' % (new_couplings[icoup],prefactors,value,all_couplings[icoup])
	elif(new_couplings[icoup]) :
	all_couplings[icoup] = new_couplings[icoup]
	# return the results
	return (ordering,all_couplings)

	def processTensorCouplings(lorentztag,vertex,model,parmsubs,all_couplings) :
	def evaluate(x):
	import cmath
	return eval(x,
	{'cmath':cmath,
	'complexconjugate':model.function_library.complexconjugate},
	parmsubs)
	# check for fermion vertices (i.e. has L/R couplings)
	fermions = "FF" in lorentztag
	# test and process the values of the couplings
	tval = [False]*3
	value = [False]*3
	# loop over the colours
	for icolor in range(0,len(all_couplings)) :
	lmax = len(all_couplings[icolor])
	if(fermions) : lmax /=3
	# loop over the different terms
	for ix in range(0,lmax) :
	test = [False]*3
	# normal case
	if( not fermions ) :
	test[0] = all_couplings[icolor][ix]
	else :
	# first case vector but no L/R couplings
	if( not all_couplings[icolor][lmax+ix] and
	not all_couplings[icolor][2*lmax+ix] ) :
	test[0] = all_couplings[icolor][ix]
	# special for mass terms and massless particles
	if(not all_couplings[icolor][ix]) :
	code = abs(vertex.particles[0].pdg_code)
	if(ix==6 and code ==12 or code ==14 or code==16) :
	continue
	else :
	raise SkipThisVertex()
	# second case L/R couplings
	elif( not all_couplings[icolor][ix] ) :
	# L=R, replace with vector
	if(all_couplings[icolor][lmax+ix] ==
	all_couplings[icolor][2*lmax+ix]) :
	test[0] = all_couplings[icolor][lmax+ix]
	else :
	test[1] = all_couplings[icolor][lmax+ix]
	test[2] = all_couplings[icolor][2*lmax+ix]
	else :
	raise SkipThisVertex()
	# special handling of mass terms
	# scalar divide by m**2
	if((ix==3 and lorentztag=="SST") or
	( ix>=10 and ix<=12 and lorentztag=="VVT" )) :
	for i in range(0,len(test)) :
	if(test[i]) :
	test[i] = '(%s)/%s**2' % (test[i],vertex.particles[0].mass.value)
	# fermions divide by 4*m
	elif(ix==6 and lorentztag=="FFT" and
	float(vertex.particles[0].mass.value) != 0. ) :
	for i in range(0,len(test)) :
	if(test[i]) :
	test[i] = '-(%s)/%s/4' % (test[i],vertex.particles[0].mass.value)
	# set values on first pass
	if(not tval[0] and not tval[1] and not tval[2]) :
	value = test
	for i in range(0,len(test)) :
	if(test[i]) : tval[i] = evaluate(test[i])
	else :
	for i in range(0,len(test)) :
	if(not test[i] and not tval[i]) :
	continue
	if(not test[i] or not tval[i]) :
	# special for mass terms and vectors
	if(lorentztag=="VVT" and ix >=10 and ix <=12 and
	float(vertex.particles[0].mass.value) == 0. ) :
	continue
	# special for vector gauge terms
	if(lorentztag=="VVT" and ix>=13) :
	continue
	raise SkipThisVertex()
	tval2 = evaluate(test[i])
	if(abs(tval[i]-tval2)>1e-6) :
	# special for fermion mass term if fermions massless
	if(lorentztag=="FFT" and ix ==6 and tval2 == 0. and
	float(vertex.particles[0].mass.value) == 0. ) :
	continue
	raise SkipThisVertex()
	# simple clean up
	for i in range(0,len(value)):
	if(value[i]) :
	value[i] = value[i].replace("(1.0) * ","").replace(" * (1)","")
	# put everything together
	- coup_left = []
	- coup_right = []
	- coup_norm = []
	+ coup_left = 0.
	+ coup_right = 0.
	+ coup_norm = 0.
	if(lorentztag == "SST" or lorentztag == "VVT" or
	lorentztag == "VVVT" or lorentztag == "FFT" ) :
	coup_norm.append(value[0])
	if(value[1] or value[2]) :
	raise SkipThisVertex()
	elif(lorentztag=="FFVT") :
	if(not value[1] and not value[2]) :
	- coup_norm.append(value[0])
	- coup_left.append("1.")
	- coup_right.append("1.")
	+ coup_norm = value[0]
	+ coup_left = "1."
	+ coup_right = "1."
	elif(not value[0]) :
	- coup_norm.append("1.")
	+ coup_norm = "1."
	if(value[1] and value[2]) :
	- coup_left.append(value[1])
	- coup_right.append(value[2])
	+ coup_left = value[1]
	+ coup_right = value[2]
	elif(value[1]) :
	- coup_left.append(value[1])
	- coup_right.append("0.")
	+ coup_left = value[1]
	+ coup_right = "0."
	elif(value[2]) :
	- coup_left.append("0.")
	- coup_right.append(value[2])
	+ coup_left = "0."
	+ coup_right = value[2]
	else :
	raise SkipThisVertex()
	else :
	raise SkipThisVertex()
	else :
	raise SkipThisVertex()
	# return the answer
	return (coup_left,coup_right,coup_norm)

	def EWVVVVCouplings(vertex,L) :
	terms=['Metric(1,2)*Metric(3,4)',
	'Metric(1,3)*Metric(2,4)',
	'Metric(1,4)*Metric(2,3)']

	structure1 = L.structure.split()

	structures =[]
	sign=''
	for struct in structure1 :
	if(struct=='+') :
	continue
	elif(struct=='-') :
	sign='-'
	else :
	structures.append(sign+struct.strip())
	sign=''
	factors=[]
	for term in terms:
	for struct in structures :
	if term in struct :
	reminder = struct.replace(term,'1.',1)
	try:
	factors.append(eval(reminder, {'cmath':cmath} ))
	except NameError:
	name_error = True
	else:
	name_error = False

	if len(factors) != 3 or name_error:
	sys.stderr.write(
	'Warning: unsupported {tag} ( {ps} ) Lorentz structure in {name}:\n{lstr}\n'
	.format(tag=unique_lorentztag(vertex), name=vertex.name,
	lstr=L.structure, ps=' '.join(map(str,vertex.particles)))
	)
	raise SkipThisVertex()

	factor=0.
	order=[]
	if(factors[0]==-2.factors[1] and factors[0]==-2.factors[2] ) :
	order=[0,1,2,3]
	factor = factors[0]/2.
	elif(factors[1]==-2.factors[0] and factors[1]==-2.factors[2] ) :
	order=[0,2,1,3]
	factor = factors[1]/2.
	elif(factors[2]==-2.factors[0] and factors[2]==-2.factors[1] ) :
	order=[0,3,1,2]
	factor = factors[2]/2.
	else:
	sys.stderr.write(
	'Warning: unsupported {tag} ( {ps} ) Lorentz structure in {name}:\n{lstr}\n'
	.format(tag=unique_lorentztag(vertex), name=vertex.name,
	lstr=L.structure, ps=' '.join(map(str,vertex.particles)))
	)
	raise SkipThisVertex()


	pattern = \
	"bool done[4]={false,false,false,false};\n" + \
	" tcPDPtr part[4]={p1,p2,p3,p4};\n" + \
	" unsigned int iorder[4]={0,0,0,0};\n" + \
	" for(unsigned int ix=0;ix<4;++ix) {\n" + \
	" if(!done[0] && part[ix]->id()==%s) {done[0]=true; iorder[%s] = ix; continue;}\n" + \
	" if(!done[1] && part[ix]->id()==%s) {done[1]=true; iorder[%s] = ix; continue;}\n" + \
	" if(!done[2] && part[ix]->id()==%s) {done[2]=true; iorder[%s] = ix; continue;}\n" + \
	" if(!done[3] && part[ix]->id()==%s) {done[3]=true; iorder[%s] = ix; continue;}\n" + \
	" }\n" + \
	" setType(2);\n" + \
	" setOrder(iorder[0],iorder[1],iorder[2],iorder[3]);"
	ordering=pattern % ( vertex.particles[0].pdg_code,order[0],
	vertex.particles[1].pdg_code,order[1],
	vertex.particles[2].pdg_code,order[2],
	vertex.particles[3].pdg_code,order[3] )
	return (ordering,factor)

	def changeSign(sign1,sign2) :
	if((sign1=="+" and sign2=="+") or\
	(sign1=="-" and sign2=="-")) :
	return "+"
	else :
	return "-"

	def epsilonOrder(eps) :
	terms = eps.strip("Epsilon(").strip(")").split(",")
	sign=1.
	for iy in range(0,len(terms)) :
	for ix in range(-1,-len(terms)+iy,-1) :
	swap = False
	if(len(terms[ix])==1 and len(terms[ix-1])==1) :
	swap = int(terms[ix])<int(terms[ix-1])
	elif(len(terms[ix])==2 and len(terms[ix-1])==2) :
	swap = int(terms[ix][1])<int(terms[ix-1][1])
	elif(len(terms[ix])==1 and len(terms[ix-1])==2) :
	swap = True
	if(swap) :
	sign *=-1.
	terms[ix],terms[ix-1] = terms[ix-1],terms[ix]
	return (sign,"Epsilon(%s,%s,%s,%s)" % (terms[0],terms[1],terms[2],terms[3]))


	def VVSEpsilon(couplings,struct) :
	if(struct.find("Epsilon")<0) :
	return
	fact=""
	sign="+"
	if(struct[-1]==")") :
	fact=struct.split("(")[0]
	if(fact.find("Epsilon")>=0) :
	fact=""
	else :
	struct=struct.split("(",1)[1][:-1]
	if(fact[0]=="-") :
	sign="-"
	fact=fact[1:]
	split = struct.split("*")
	# find the epsilon
	eps=""
	for piece in split :
	if(piece.find("Epsilon")>=0) :
	eps=piece
	split.remove(piece)
	break
	# and any prefactors
	for piece in split :
	if(piece.find("P(")<0) :
	split.remove(piece)
	if(piece[0]=="+" or piece[0]=="-") :
	sign=changeSign(sign,piece[0])
	piece=piece[1:]
	if(fact=="") :
	fact=piece
	else :
	fact = "( %s ) * ( %s )" % (fact , piece)
	# now sort out the momenta
	for piece in split :
	terms=piece.split(",")
	terms[0]=terms[0].strip("P(")
	terms[1]=terms[1].strip(")")
	eps=eps.replace(terms[0],"P%s"%terms[1])
	(nsign,eps)=epsilonOrder(eps)
	if(nsign<0) : sign=changeSign(sign,"-")
	if(fact=="") : fact="1."
	if(eps!="Epsilon(1,2,P1,P2)") :
	return
	if(couplings[6]==0.) :
	couplings[6] = "( %s%s )" % (sign,fact)
	else :
	couplings[6] = "( %s ) + ( %s%s )" % (couplings[6],sign,fact)

	-def VVSCouplings(vertex,coupling,prefactors,L,lorentztag) :
	- # split the structure into its different terms for analysis
	- ordering=""
	+
	+def scalarVectorCouplings(vertex,value,prefactors,L,lorentztag,pos,
	+ all_couplings,append,kinematics) :
	+ # set up the types of term we are looking for
	+ if(lorentztag=="VVS") :
	+ couplings=[0.,0.,0.,0.,0.,0.,0.]
	+ terms=[['P(-1,1)','P(-1,2)','Metric(1,2)'],
	+ ['P(1,1)','P(2,1)'],
	+ ['P(1,1)','P(2,2)'],
	+ ['P(1,2)','P(2,1)'],
	+ ['P(1,2)','P(2,2)'],
	+ ['Metric(1,2)']]
	+ elif(lorentztag=="VVSS") :
	+ couplings=[0.]
	+ terms=[['Metric(1,2)']]
	+ elif(lorentztag=="VSS"):
	+ couplings=[0.,0.]
	+ terms=[['P(1,3)'],['P(1,2)']]
	+ # extract the lorentz structures
	structure1 = L.structure.split()
	structures =[]
	sign=''
	- # extract the lorentz structures
	for struct in structure1 :
	if(struct=='+') :
	continue
	elif(struct=='-') :
	sign='-'
	else :
	structures.append(sign+struct.strip())
	sign=''
	# handle the scalar couplings
	- couplings=[0.,0.,0.,0.,0.,0.,0.]
	- terms=[['P(-1,1)','P(-1,2)','Metric(1,2)'],
	- ['P(1,1)','P(2,1)'],
	- ['P(1,1)','P(2,2)'],
	- ['P(1,2)','P(2,1)'],
	- ['P(1,2)','P(2,2)'],['Metric(1,2)']]
	itype=-1
	for term in terms:
	itype+=1
	for perm in itertools.permutations(term):
	label = '*'.join(perm)
	for istruct in range(0,len(structures)) :
	if label in structures[istruct] :
	reminder = structures[istruct].replace(label,'1.',1)
	couplings[itype]+=eval(reminder, {'cmath':cmath} )
	structures[istruct]='Done'
	+ # special for VVS and epsilon
	# handle the pseudoscalar couplings
	for struct in structures :
	if(struct != "Done" ) :
	- VVSEpsilon(couplings,struct)
	- # evaluate the prefactor
	- if type(coupling) is not list:
	- value = coupling.value
	- else:
	- value = "("
	- for coup in coupling :
	- value += '+(%s)' % coup.value
	- value +=")"
	+ if(lorentztag=="VVS") :
	+ VVSEpsilon(couplings,struct)
	+ else :
	+ raise SkipThisVertex()
	# put it all together
	- for ic in range(0,len(couplings)) :
	- if(couplings[ic]!=0.) :
	- couplings[ic] = '(%s) * (%s) * (%s)' % (prefactors,value,couplings[ic])
	- return couplings
	+ if(len(all_couplings)==0) :
	+ for ic in range(0,len(couplings)) :
	+ if(couplings[ic]!=0.) :
	+ all_couplings.append('(%s) * (%s) * (%s)' % (prefactors,value,couplings[ic]))
	+ else :
	+ all_couplings.append(False)
	+ else :
	+ for ic in range(0,len(couplings)) :
	+ if(couplings[ic]!=0. and all_couplings[ic]) :
	+ all_couplings[ic] = '(%s) * (%s) * (%s) + (%s) ' % (prefactors,value,
	+ couplings[ic],all_couplings[ic])
	+ elif(couplings[ic]!=0) :
	+ all_couplings[ic] = '(%s) * (%s) * (%s) ' % (prefactors,value,couplings[ic])
	+ return all_couplings
	+
	+def processScalarVectorCouplings(lorentztag,vertex,model,parmsubs,all_couplings,kinematics) :
	+ def evaluate(x):
	+ import cmath
	+ return eval(x,
	+ {'cmath':cmath,
	+ 'complexconjugate':model.function_library.complexconjugate},
	+ parmsubs)
	+ # check the values
	+ tval = [False]*len(all_couplings[0])
	+ value =[False]*len(all_couplings[0])
	+ for icolor in range(0,len(all_couplings)) :
	+ for ix in range(0,len(all_couplings[icolor])) :
	+ if(not value[ix]) :
	+ value[ix] = all_couplings[icolor][ix]
	+ if(value[ix] and not tval[ix]) :
	+ tval[ix] = evaluate(value[ix])
	+ elif(value[ix]) :
	+ tval2 = evaluate(all_couplings[icolor][0])
	+ if(abs(tval[ix]-tval2)>1e-6) :
	+ raise SkipThisVertex()
	+
	+ append = ""
	+ symbols = set()
	+ coup_norm=0.
	+ if(lorentztag=="VVS") :
	+ if(not value[0] and not value[1] and not value[2] and \
	+ not value[3] and not value[4] and not value[6] and value[5]) :
	+ coup_norm=value[5]
	+ else :
	+ for ix in range(0,len(value)) :
	+ if(value[ix]) :
	+ value[ix], sym = py2cpp(value[ix])
	+ symbols \|= sym
	+ else :
	+ value[ix]=0.
	+ lorentztag = 'GeneralVVS'
	+ kinematics='true'
	+ # g_mu,nv piece of coupling
	+ if(value[5]!=0.) :
	+ append +='a00( %s + Complex(( %s )* GeV2/invariant(1,2)));\n' % ( value[0],value[5])
	+ else :
	+ append +='a00( %s );\n' % value[0]
	+ # other couplings
	+ append += 'a11( %s );\n a12( %s );\n a21( %s );\n a22( %s );\n aEp( %s );\n' % \
	+ ( value[1],value[2],value[3],value[4],value[6] )
	+ coup_norm="1."
	+ elif(lorentztag=="VVSS") :
	+ coup_norm = value[0]
	+ elif(lorentztag=="VSS") :
	+ if(abs(tval[0]+tval[1])>1e-6) :
	+ raise SkipThisVertex()
	+ coup_norm = value[1]
	+ append = 'if(p2->id()!=%s){norm(-norm());}' \
	+ % vertex.particles[1].pdg_code
	+ # # cleanup and return the answer
	+ # value = value.replace("(1.0) * ","").replace(" * (1)","")
	+ # return [value]
	+ return (coup_norm,append,lorentztag,kinematics,symbols)

	def getIndices(term) :
	if(term[0:2]=="P(") :
	indices = term.strip(")").strip("P(").split(",")
	mom = int(indices[1])
	index = int(indices[0])
	return (True,mom,index)
	else :
	return (False,0,0)


	def lorentzScalar(vertex,L) :
	dotProduct = """\
	0.5*(invariant( i[{i1}], i[{i2}] ) - invariant( i[{i1}], i[{i1}] ) -invariant( i[{i2}] , i[{i2}] ))/GeV2"""
	structure = L.structure.split("*")
	worked = False
	mom=-1
	output=""
	while True :
	term = structure[-1]
	structure.pop()
	(momentum,mom,index) = getIndices(term)
	if( not momentum) : break
	# look for the matching momenta
	for term in structure :
	(momentum,mom2,index2) = getIndices(term)
	if(index2==index) :
	structure.remove(term)
	dot = dotProduct.format(i1=mom-1,i2=mom2-1)
	if(output=="") :
	output = dot
	else :
	output = " ( %s) * ( %s ) " (output,dot)
	if(len(structure)==0) :
	worked = True
	break
	if( not worked ) :
	return False
	else :
	return output

	kinematicsline = """\
	long id [3]={{{id1},{id2},{id3}}};
	long id2[3]={{p1->id(),p2->id(),p3->id()}};
	unsigned int i[3];
	for(unsigned int ix=0;ix<3;++ix) {{
	for(unsigned int iy=0;iy<3;++iy) {{
	if(id[ix]==id2[iy]) {{
	i[ix] = iy;
	id2[iy]=0;
	break;
	}}
	}}
	}}
	double hw_kine1 = {kine};
	"""

	kinematicsline2 = """\
	long id [4]={{{id1},{id2},{id3},{id4}}};
	long id2[4]={{p1->id(),p2->id(),p3->id(),p4->id()}};
	unsigned int i[4];
	for(unsigned int ix=0;ix<4;++ix) {{
	for(unsigned int iy=0;iy<4;++iy) {{
	if(id[ix]==id2[iy]) {{
	i[ix] = iy;
	id2[iy]=0;
	break;
	}}
	}}
	}}
	double hw_kine1 = {kine};
	"""

	kinematicsline3 ="""\
	double hw_kine{i} = {kine};
	"""

	def scalarCouplings(vertex,value,prefactors,L,lorentztag,pos,
	all_couplings,ordering,kinematics) :
	try :
	val = int(L.structure)
	except :
	output = lorentzScalar(vertex,L)
	if( not output ) :
	raise SkipThisVertex()
	else :
	if(ordering=="") :
	if(lorentztag=="SSS") :
	ordering = kinematicsline.format(id1=vertex.particles[0].pdg_code,
	id2=vertex.particles[1].pdg_code,
	id3=vertex.particles[2].pdg_code,
	kine=output)
	else :
	ordering = kinematicsline2.format(id1=vertex.particles[0].pdg_code,
	id2=vertex.particles[1].pdg_code,
	id3=vertex.particles[2].pdg_code,
	id4=vertex.particles[2].pdg_code,
	kine=output)
	value = "(%s) *(hw_kine1)" % value
	else :
	osplit=ordering.split("\n")
	i=-1
	while osplit[i]=="":
	i=i-1
	ikin=int(osplit[i].split("=")[0].replace("double hw_kine",""))+1
	ordering +=kinematicsline3.format(kine=output,i=ikin)
	value = "(%s) *(hw_kine%s)" % (value,ikin)
	kinematics="true"
	if(len(all_couplings)==0) :
	all_couplings.append('(%s) * (%s)' % (prefactors,value))
	else :
	all_couplings[0] = '(%s) * (%s) + (%s)' % (prefactors,value)
	return (ordering, kinematics,all_couplings)

	def processScalarCouplings(lorentztag,vertex,model,parmsubs,all_couplings) :
	def evaluate(x):
	import cmath
	return eval(x,
	{'cmath':cmath,
	'complexconjugate':model.function_library.complexconjugate},
	parmsubs)
	tval = False
	value = False
	for icolor in range(0,len(all_couplings)) :
	if(len(all_couplings[icolor])!=1) :
	raise SkipThisVertex
	if(not value) :
	value = all_couplings[icolor][0]
	m = re.findall('hw_kine[0-9]*', all_couplings[icolor][0])
	if m:
	for kine in m:
	# bizarre number for checks, must be a better option
	parmsubs[kine] = 987654321.
	if(not tval) :
	tval = evaluate(value)
	else :
	tval2 = evaluate(all_couplings[icolor][0])
	if(abs(tval[i]-tval2)>1e-6) :
	raise SkipThisVertex()
	# cleanup and return the answer
	- value = value.replace("(1.0) * ","").replace(" * (1)","")
	- return [value]
	+ return value.replace("(1.0) * ","").replace(" * (1)","")
	diff --git a/Models/Feynrules/python/ufo2peg/vertices.py b/Models/Feynrules/python/ufo2peg/vertices.py
	--- a/Models/Feynrules/python/ufo2peg/vertices.py
	+++ b/Models/Feynrules/python/ufo2peg/vertices.py
	@@ -1,627 +1,610 @@
	import sys,pprint

	from .helpers import CheckUnique,getTemplate,writeFile,qcd_qed_orders,def_from_model
	from .lorentzparser import parse_lorentz
	from .converter import py2cpp
	from .collapse_vertices import collapse_vertices
	-from .check_lorentz import tensorCouplings,VVVordering,VVSCouplings,EWVVVVCouplings,lorentzScalar,\
	- processTensorCouplings,scalarCouplings,processScalarCouplings
	+from .check_lorentz import tensorCouplings,VVVordering,EWVVVVCouplings,lorentzScalar,\
	+ processTensorCouplings,scalarCouplings,processScalarCouplings,scalarVectorCouplings,\
	+ processScalarVectorCouplings
	from .helpers import SkipThisVertex

	# names of goldstone bosons
	gsnames = ['goldstone','goldstoneboson','GoldstoneBoson']

	# prefactors for vertices
	lfactors = {
	'FFV' : '-complex(0,1)', # ok
	'VVV' : 'complex(0,1)', # changed to fix ttbar
	'VVVS' : 'complex(0,1)', # should be as VVV
	'VVVV' : 'complex(0,1)',
	'VVS' : '-complex(0,1)',
	'VSS' : '-complex(0,1)', # changed to minus to fix dL ->x1 W- d
	'SSS' : '-complex(0,1)', # ok
	'VVSS' : '-complex(0,1)', # ok
	'VVT' : 'complex(0,2)',
	'VVVT' : '-complex(0,2)',
	'SSSS' : '-complex(0,1)', # ok
	'FFS' : '-complex(0,1)', # ok
	'SST' : 'complex(0,2)',
	'FFT' : '-complex(0,8)',
	'FFVT' : '-complex(0,4)',
	}

	# template for the header for a vertex
	VERTEXHEADER = """\
	#include "ThePEG/Helicity/Vertex/{spindirectory}/{lorentztag}Vertex.h"
	"""

	# template for the implmentation for a vertex
	VERTEXCLASS = getTemplate('Vertex_class')

	# template for the .cc file for vertices
	VERTEX = getTemplate('Vertex.cc')

	vertexline = """\
	create Herwig::{modelname}V_{vname} /Herwig/{modelname}/V_{vname}
	insert {modelname}:ExtraVertices 0 /Herwig/{modelname}/V_{vname}
	"""


	def get_lorentztag(spin):
	"""Produce a ThePEG spin tag for the given numeric FR spins."""
	spins = { 1 : 'S', 2 : 'F', 3 : 'V', -1 : 'U', 5 : 'T' }
	result = [ spins[s] for s in spin ]

	def spinsort(a,b):
	"""Helper function for ThePEG's FVST spin tag ordering."""
	if a == b: return 0
	for letter in 'UFVST':
	if a == letter: return -1
	if b == letter: return 1

	result = sorted(result, cmp=spinsort)
	return ''.join(result)

	def unique_lorentztag(vertex):
	"""Check and return the Lorentz tag of the vertex."""
	unique = CheckUnique()
	for l in vertex.lorentz:
	lorentztag = get_lorentztag(l.spins)
	unique( lorentztag )
	lname = l.name[:len(lorentztag)]
	if sorted(lorentztag) != sorted(lname):
	raise Exception("Lorentztags: %s is not %s in %s"
	% (lorentztag,lname,vertex))
	# if lorentztag != lname:
	# sys.stderr.write("Warning: Lorentz tag ordering: %s is not %s in %s\n"
	# % (lorentztag,lname,vertex))

	return lorentztag

	def colors(vertex) :
	try:
	unique = CheckUnique()
	for pl in vertex.particles_list:
	struct = [ p.color for p in pl ]
	unique(struct)
	except:
	struct = [ p.color for p in vertex.particles ]
	pos = colorpositions(struct)
	L = len(struct)
	return (L,pos)

	def colorfactor(vertex,L,pos):
	def match(patterns):
	result = [ p == t
	for p,t in zip(patterns,vertex.color) ]
	return all(result)

	label = None
	l = lambda c: len(pos[c])
	if l(1) == L:
	label = ('1',)
	if match(label): return ('1',)

	elif l(3) == l(-3) == 1 and l(1) == L-2:
	nums = [pos[3][0], pos[-3][0]]
	label = ('Identity({0},{1})'.format(*sorted(nums)),)
	if match(label): return ('1',)

	elif l(6) == l(-6) == 1 and l(1) == L-2:
	nums = [pos[6][0], pos[-6][0]]
	label = ('Identity({0},{1})'.format(*sorted(nums)),)
	if match(label): return ('1',)

	elif l(6) == l(-6) == 2 and L==4:
	sys.stderr.write(
	'Warning: Unknown colour structure 6 6 6~ 6~ ( {ps} ) in {name}.\n'
	.format(name=vertex.name, ps=' '.join(map(str,vertex.particles)))
	)
	raise SkipThisVertex()

	elif l(8) == l(3) == l(-3) == 1 and l(1) == L-3:
	label = ('T({g},{q},{qb})'.format(g=pos[8][0],q=pos[3][0],qb=pos[-3][0]),)
	if match(label): return ('1',)

	elif l(8) == l(6) == l(-6) == 1 and l(1) == L-3:
	label = ('T6({g},{s},{sb})'.format(g=pos[8][0],s=pos[6][0],sb=pos[-6][0]),)
	if match(label): return ('1',)

	elif l(6) == 1 and l(-3) == 2 and L==3:
	label = ('K6({s},{qb1},{qb2})'.format(s=pos[6][0],qb1=pos[-3][0],qb2=pos[-3][1]),)
	if match(label): return ('1',)

	elif l(-6) == 1 and l(3) == 2 and L==3:
	label = ('K6Bar({sb},{q1},{q2})'.format(sb=pos[-6][0],q1=pos[3][0],q2=pos[3][1]),)
	if match(label): return ('1',)

	elif l(3) == L == 3:
	label = ('Epsilon(1,2,3)',)
	if match(label): return ('1',) # TODO check factor!

	elif l(-3) == L == 3:
	label = ('EpsilonBar(1,2,3)',)
	if match(label): return ('1',) # TODO check factor!

	elif l(8) == L == 3:
	# if lorentz is FFV get extra minus sign
	lorentztag = unique_lorentztag(vertex)
	factor = '*(-1)' if lorentztag in ['FFV'] else ''
	label = ('f(1,2,3)',)
	if match(label): return ('-complex(0,1)%s'%factor,)
	label = ('f(3,2,1)',)
	if match(label): return ('complex(0,1)%s'%factor,)
	label = ('f(2,1,3)',)
	if match(label): return ('complex(0,1)%s'%factor,)

	elif l(8) == L == 4:
	label = ('f(-1,1,2)*f(3,4,-1)',
	'f(-1,1,3)*f(2,4,-1)',
	'f(-1,1,4)*f(2,3,-1)',
	)
	if match(label): return ('-1./3.','-1./3.','-1./3.')

	elif l(8) == 2 and l(3) == l(-3) == 1 and L==4:
	subs = {
	'g1' : pos[8][0],
	'g2' : pos[8][1],
	'qq' : pos[3][0],
	'qb' : pos[-3][0]
	}
	label = ('T({g1},-1,{qb})T({g2},{qq},-1)'.format(*subs),
	'T({g1},{qq},-1)T({g2},-1,{qb})'.format(*subs))
	if match(label): return ('0.5','0.5')

	elif l(8) == 2 and l(6) == l(-6) == 1 and L==4:
	subs = {
	'g1' : pos[8][0],
	'g2' : pos[8][1],
	'qq' : pos[6][0],
	'qb' : pos[-6][0]
	}
	label = ('T6({g1},-1,{qb})T6({g2},{qq},-1)'.format(*subs),
	'T6({g1},{qq},-1)T6({g2},-1,{qb})'.format(*subs))
	if match(label): return ('0.5','0.5')

	elif l(8) == 2 and l(8)+l(1)==L :
	subs = { 'g1' : pos[8][0], 'g2' : pos[8][1] }
	label = ('Identity({g1},{g2})'.format(**subs),)
	if match(label) : return ('1.',)

	elif l(8) == 3 and l(1)==1 and L==4 :
	label = ('f(1,2,3)',)
	if match(label): return ('-complex(0,1)',)
	label = ('f(3,2,1)',)
	if match(label): return ('complex(0,1)',)
	label = ('f(2,1,3)',)
	if match(label): return ('complex(0,1)',)

	sys.stderr.write(
	"Warning: Unknown colour structure {color} ( {ps} ) in {name}.\n"
	.format(color = ' '.join(vertex.color), name = vertex.name,
	ps = ' '.join(map(str,vertex.particles)))
	)
	raise SkipThisVertex()

	def colorpositions(struct):
	positions = {
	1 : [],
	3 : [],
	-3 : [],
	6 : [],
	-6 : [],
	8 : [],
	}
	for i,s in enumerate(struct,1):
	positions[s].append(i)
	return positions

	def spindirectory(lt):
	"""Return the spin directory name for a given Lorentz tag."""
	if 'T' in lt:
	spin_directory = 'Tensor'
	elif 'S' in lt:
	spin_directory = 'Scalar'
	elif 'V' in lt:
	spin_directory = 'Vector'
	else:
	raise Exception("Unknown Lorentz tag {lt}.".format(lt=lt))
	return spin_directory

	def write_vertex_file(subs):
	'Write the .cc file for some vertices'
	newname = '%s_Vertices_%03d.cc' % (subs['ModelName'],subs['vertexnumber'])
	subs['newname'] = newname
	writeFile( newname, VERTEX.substitute(subs) )

	def checkGhostGoldstoneVertex(lorentztag,vertex) :
	'check if vertex has ghosts or goldstones'
	# remove vertices involving ghost fields
	if 'U' in lorentztag:
	return True
	# remove vertices involving goldstones
	for p in vertex.particles:
	def gstest(name):
	try:
	return getattr(p,name)
	except AttributeError:
	return False
	if any(map(gstest, gsnames)):
	return True
	return False

	def calculatePrefactor(globalsign,lorentztag,lf,cf) :
	if(globalsign!=1.) :
	prefactors = '(%s) * (%s) * (%s)' \
	% (globalsign**(len(lorentztag)-2),lf,cf)
	else :
	prefactors = '(%s) * (%s)' \
	% (lf,cf)
	return prefactors
	# fact=[]
	# if(globalsign!=1.) :
	# fact.append(globalsign**(len(lorentztag)-2))
	# if(lf!="1") :
	# fact.append(lf)
	# if(cf!="1") :
	# fact.append(cf)
	# if(len(fact)==0) : return "1"
	# prefactor = '(%s)' % fact[0]
	# for ix in range(1,len(fact)) :
	# prefactor = '%s * (%s)' % (prefactor,fact[ix])
	# return prefactor

	def couplingValue(coupling) :
	if type(coupling) is not list:
	value = coupling.value
	else:
	value = "("
	for coup in coupling :
	value += '+(%s)' % coup.value
	value +=")"
	return value

	class VertexConverter:
	'Convert the vertices in a FR model to extract the information ThePEG needs.'
	def __init__(self,model,parmsubs) :
	'Initialize the parameters'
	self.ONE_EACH=True
	self.verbose=False
	self.vertex_skipped=False
	self.ignore_skipped=False
	self.model=model
	self.all_vertices= []
	self.modelname=""
	self.globalsign=self.global_sign()
	self.no_generic_loop_vertices = False
	self.parmsubs = parmsubs

	def global_sign(self):
	'Initial pass to find global sign at the moment does nothing'
	return 1.0
	# for v in self.model.all_vertices:
	# pids = sorted([ p.pdg_code for p in v.particles ])
	# if pids != [-11,11,22]: continue
	# coup = v.couplings
	# assert( len(coup) == 1 )
	# val = coup.values()[0].value
	# val = evaluate(val)
	# assert( val.real == 0 )
	# return 1 if val.imag > 0 else -1

	def readArgs(self,args) :
	'Extract the relevant command line arguments'
	self.ignore_skipped = args.ignore_skipped
	self.verbose = args.verbose
	self.modelname = args.name
	self.no_generic_loop_vertices = args.no_generic_loop_vertices

	def should_print(self) :
	'Check if we should output the results'
	return not self.vertex_skipped or self.ignore_skipped

	def convert(self) :
	'Convert the vertices'
	if(self.verbose) :
	print 'verbose mode on: printing all vertices'
	print '-'*60
	labels = ('vertex', 'particles', 'Lorentz', 'C_L', 'C_R', 'norm')
	pprint.pprint(labels)
	# check if we should merge vertices
	if(self.ONE_EACH) :
	self.all_vertices = self.model.all_vertices
	else:
	self.all_vertices = collapse_vertices(self.model.all_vertices)
	# convert the vertices
	vertexclasses, vertexheaders = [], set()
	for vertexnumber,vertex in enumerate(self.all_vertices,1) :
	# process the vertex
	(skip,vertexClass,vertexHeader) = \
	self.processVertex(vertexnumber,vertex)
	# check it can be handled
	if(skip) : continue
	# add to the list
	vertexclasses.append(vertexClass)
	vertexheaders.add(vertexHeader)
	WRAP = 25
	if vertexnumber % WRAP == 0:
	write_vertex_file({'vertexnumber' : vertexnumber//WRAP,
	'vertexclasses' : '\n'.join(vertexclasses),
	'vertexheaders' : ''.join(vertexheaders),
	'ModelName' : self.modelname})
	vertexclasses = []
	vertexheaders = set()
	# exit if there's vertices we can't handle
	if not self.should_print():
	sys.stderr.write(
	"""
	Error: The conversion was unsuccessful, some vertices could not be
	generated. If you think the missing vertices are not important
	and want to go ahead anyway, use --ignore-skipped.
	Herwig may not give correct results, though.
	"""
	)
	sys.exit(1)
	# if still stuff to output it do it
	if vertexclasses:
	write_vertex_file({'vertexnumber' : vertexnumber//WRAP + 1,
	'vertexclasses' : '\n'.join(vertexclasses),
	'vertexheaders' : ''.join(vertexheaders),
	'ModelName' : self.modelname})

	print '='*60

	+ def setCouplingPtrs(self,lorentztag,qcd,append,prepend) :
	+ couplingptrs = [',tcPDPtr']*len(lorentztag)
	+ if lorentztag == 'VSS':
	+ couplingptrs[1] += ' p2'
	+ elif lorentztag == 'FFV':
	+ couplingptrs[0] += ' p1'
	+ elif (lorentztag == 'VVV' or lorentztag == 'VVVS' or
	+ lorentztag == "SSS" or lorentztag == "VVVT" ) \
	+ and (append or prepend ) :
	+ couplingptrs[0] += ' p1'
	+ couplingptrs[1] += ' p2'
	+ couplingptrs[2] += ' p3'
	+ elif (lorentztag == 'VVVV' and qcd != 2) or\
	+ (lorentztag == "SSSS" and prepend ):
	+ couplingptrs[0] += ' p1'
	+ couplingptrs[1] += ' p2'
	+ couplingptrs[2] += ' p3'
	+ couplingptrs[3] += ' p4'
	+ return couplingptrs
	+
	def processVertex(self,vertexnumber,vertex) :
	# get the Lorentz tag for the vertex
	lorentztag = unique_lorentztag(vertex)
	# check if we should skip the vertex
	vertex.herwig_skip_vertex = checkGhostGoldstoneVertex(lorentztag,vertex)
	if(vertex.herwig_skip_vertex) :
	return (True,"","")
	# get the factor for the vertex
	try:
	lf = lfactors[lorentztag]
	except KeyError:
	msg = 'Warning: Lorentz structure {tag} ( {ps} ) in {name} ' \
	'is not supported.\n'.format(tag=lorentztag, name=vertex.name,
	ps=' '.join(map(str,vertex.particles)))
	sys.stderr.write(msg)
	vertex.herwig_skip_vertex = True
	self.vertex_skipped=True
	return (True,"","")
	# get the ids of the particles at the vertex
	if self.ONE_EACH:
	plistarray = [ ','.join([ str(p.pdg_code) for p in vertex.particles ]) ]
	else:
	plistarray = [ ','.join([ str(p.pdg_code) for p in pl ])
	for pl in vertex.particles_list ]
	# parse the colour structure for the vertex
	try:
	L,pos = colors(vertex)
	cf = colorfactor(vertex,L,pos)
	except SkipThisVertex:
	msg = 'Warning: Color structure for vertex ( {ps} ) in {name} ' \
	'is not supported.\n'.format(tag=lorentztag, name=vertex.name,
	ps=' '.join(map(str,vertex.particles)))
	sys.stderr.write(msg)
	vertex.herwig_skip_vertex = True
	self.vertex_skipped=True
	return (True,"","")

	### classname
	classname = 'V_%s' % vertex.name
	# try to extract the couplings
	try:
	- (coup_left,coup_right,coup_norm,couplings_VVS,kinematics,qcd,qed,prepend,append) = \
	+ (all_couplings,kinematics,qcd,qed,prepend,append) = \
	self.extractCouplings(lorentztag,pos,lf,cf,vertex)
	except SkipThisVertex:
	msg = 'Warning: Lorentz structure {tag} ( {ps} ) in {name} ' \
	'is not supported, may have a non-perturbative form.\n'.format(tag=lorentztag, name=vertex.name,
	ps=' '.join(map(str,vertex.particles)))
	sys.stderr.write(msg)
	vertex.herwig_skip_vertex = True
	self.vertex_skipped=True
	return (True,"","")

	- leftcontent = ' + '.join(coup_left) if coup_left else '0'
	- rightcontent = ' + '.join(coup_right) if coup_right else '0'
	- normcontent = ' + '.join(coup_norm) if coup_norm else '1'
	+ # set the coupling ptrs in the setCoupling call
	+ couplingptrs = self.setCouplingPtrs(lorentztag,qcd,append != '',prepend != '')

	- #print 'Left:',leftcontent
	- #print 'Right:',rightcontent
	- #print 'Norm:',normcontent
	- #print '---------------'
	- couplingptrs = [',tcPDPtr']*len(lorentztag)
	- if lorentztag == 'VSS':
	- couplingptrs[1] += ' p2'
	- elif lorentztag == 'FFV':
	- couplingptrs[0] += ' p1'
	- elif (lorentztag == 'VVV' or lorentztag == 'VVVS' or
	- lorentztag == "SSS" or lorentztag == "VVVT" ) \
	- and (append != '' or prepend != '') :
	- couplingptrs[0] += ' p1'
	- couplingptrs[1] += ' p2'
	- couplingptrs[2] += ' p3'
	- elif (lorentztag == 'VVVV' and qcd != 2) or\
	- (lorentztag == "SSSS" and prepend !=''):
	- couplingptrs[0] += ' p1'
	- couplingptrs[1] += ' p2'
	- couplingptrs[2] += ' p3'
	- couplingptrs[3] += ' p4'
	+ # final processing of the couplings
	+ symbols = set()
	+ if('T' in lorentztag) :
	+ (leftcontent,rightcontent,normcontent) = processTensorCouplings(lorentztag,vertex,self.model,
	+ self.parmsubs,all_couplings)
	+ elif(lorentztag=="SSS" or lorentztag=="SSSS") :
	+ normcontent = processScalarCouplings(lorentztag,vertex,self.model,self.parmsubs,all_couplings)
	+ elif(lorentztag=="VVS" or lorentztag =="VVSS" or lorentztag=="VSS") :
	+ normcontent,append,lorentztag,kinematics,sym = \
	+ processScalarVectorCouplings(lorentztag,vertex,self.model,self.parmsubs,all_couplings,kinematics)
	+ symbols \|=sym
	+ else :
	+ normcontent = all_couplings[0]
	+ leftcontent = all_couplings[1]
	+ rightcontent = all_couplings[2]
	+
	+

	### do we need left/right?
	- symbols = set()
	if 'FF' in lorentztag and lorentztag != "FFT":
	#leftcalc = aStoStrongCoup(py2cpp(leftcontent)[0], paramstoreplace_, paramstoreplace_expressions_)
	#rightcalc = aStoStrongCoup(py2cpp(rightcontent)[0], paramstoreplace_, paramstoreplace_expressions_)
	leftcalc, sym = py2cpp(leftcontent)
	symbols \|= sym
	rightcalc, sym = py2cpp(rightcontent)
	symbols \|= sym
	left = 'left(' + leftcalc + ');'
	right = 'right(' + rightcalc + ');'
	- elif lorentztag == 'VVS':
	- if(couplings_VVS[0]==0. and couplings_VVS[1]==0. and couplings_VVS[2]==0. and \
	- couplings_VVS[3]==0. and couplings_VVS[4]==0. and couplings_VVS[5]!=0) :
	- normcontent = couplings_VVS[5]
	- left=''
	- right=''
	- else :
	- for ix in range(0,len(couplings_VVS)) :
	- if(couplings_VVS[ix] !=0.) :
	- couplings_VVS[ix], sym = py2cpp(couplings_VVS[ix])
	- symbols \|= sym
	- lorentztag = 'GeneralVVS'
	- kinematics='true'
	- # g_mu,nv piece of coupling
	- if(couplings_VVS[5]!=0.) :
	- left = 'a00( %s + Complex(( %s )* GeV2/invariant(1,2)));' % ( couplings_VVS[0],couplings_VVS[5])
	- else :
	- left = 'a00( %s );' % couplings_VVS[0]
	- # other couplings
	- right = 'a11( %s );\n a12( %s );\n a21( %s );\n a22( %s ); aEp( %s ); ' % \
	- ( couplings_VVS[1],couplings_VVS[2],couplings_VVS[3],couplings_VVS[4],couplings_VVS[6] )
	else:
	left = ''
	right = ''
	leftcalc = ''
	rightcalc = ''
	#normcalc = aStoStrongCoup(py2cpp(normcontent)[0], paramstoreplace_, paramstoreplace_expressions_)
	normcalc, sym = py2cpp(normcontent)
	symbols \|= sym
	# UFO is GeV by default (?)
	if lorentztag in ['VVS','SSS']:
	normcalc = 'Complex((%s) * GeV / UnitRemoval::E)' % normcalc
	elif lorentztag in ['GeneralVVS']:
	normcalc = 'Complex(-(%s) * UnitRemoval::E / GeV )' % normcalc
	elif lorentztag in ['FFT','VVT', 'SST', 'FFVT', 'VVVT' , 'VVVS' ]:
	normcalc = 'Complex((%s) / GeV * UnitRemoval::E)' % normcalc
	norm = 'norm(' + normcalc + ');'
	-
	# define unkown symbols from the model
	symboldefs = [ def_from_model(self.model,s) for s in symbols ]
	### assemble dictionary and fill template
	subs = { 'lorentztag' : lorentztag, # ok
	'classname' : classname, # ok
	'symbolrefs' : '\n '.join(symboldefs),
	'left' : left, # doesn't always exist in base
	'right' : right, # doesn't always exist in base
	'norm' : norm, # needs norm, too

	#################### need survey which different setter methods exist in base classes

	'addToPlist' : '\n'.join([ 'addToList(%s);'%s for s in plistarray]),
	'parameters' : '',
	'setCouplings' : '',
	'qedorder' : qed,
	'qcdorder' : qcd,
	'couplingptrs' : ''.join(couplingptrs),
	'spindirectory' : spindirectory(lorentztag),
	'ModelName' : self.modelname,
	'prepend' : prepend,
	'append' : append,
	'kinematics' : kinematics
	} # ok

	# print info if required
	if self.verbose:
	print '-'*60
	pprint.pprint(( classname, plistarray, leftcalc, rightcalc, normcalc ))

	return (False,VERTEXCLASS.substitute(subs),VERTEXHEADER.format(**subs))

	def get_vertices(self,libname):
	vlist = ['library %s\n' % libname]
	for v in self.all_vertices:
	if v.herwig_skip_vertex: continue
	vlist.append( vertexline.format(modelname=self.modelname, vname=v.name) )
	if( not self.no_generic_loop_vertices) :
	vlist.append('insert {modelname}:ExtraVertices 0 /Herwig/{modelname}/V_GenericHPP\n'.format(modelname=self.modelname) )
	vlist.append('insert {modelname}:ExtraVertices 0 /Herwig/{modelname}/V_GenericHGG\n'.format(modelname=self.modelname) )
	return ''.join(vlist)


	def extractCouplings(self,lorentztag,pos,lf,cf,vertex) :
	coup_left = []
	coup_right = []
	coup_norm = []
	- couplings_VVS = []
	kinematics = "false"
	qcd=0
	qed=0
	prepend=""
	append=""
	unique_qcd = CheckUnique()
	unique_qed = CheckUnique()
	maxColour=0
	for (color_idx,lorentz_idx),coupling in vertex.couplings.iteritems():
	maxColour=max(maxColour,color_idx)
	all_couplings=[]
	for ix in range(0,maxColour+1) :
	all_couplings.append([])
	for colour in range(0,maxColour+1) :
	for (color_idx,lorentz_idx),coupling in vertex.couplings.iteritems() :
	if(color_idx!=colour) : continue
	qcd, qed = qcd_qed_orders(vertex, coupling)
	unique_qcd( qcd )
	unique_qed( qed )
	L = vertex.lorentz[lorentz_idx]
	prefactors = calculatePrefactor(self.globalsign,lorentztag,lf,cf[color_idx])
	# calculate the value of the coupling
	value = couplingValue(coupling)
	# handling of the different types of couplings
	if lorentztag in ['FFS','FFV']:
	left,right = parse_lorentz(L.structure)
	if left:
	coup_left.append('(%s) * (%s) * (%s)' % (prefactors,left,value))
	if right:
	coup_right.append('(%s) * (%s) * (%s)' % (prefactors,right,value))
	if lorentztag == 'FFV':
	append = ('if(p1->id()!=%s) {Complex ltemp=left(), rtemp=right(); left(-rtemp); right(-ltemp);}'
	% vertex.particles[0].pdg_code)
	elif 'T' in lorentztag :
	append, all_couplings[color_idx] = tensorCouplings(vertex,value,prefactors,L,lorentztag,pos,
	all_couplings[color_idx])
	- elif lorentztag == 'VVS' :
	- tc = VVSCouplings(vertex,coupling,prefactors,L,lorentztag)
	- if(len(couplings_VVS)==0) :
	- couplings_VVS=tc
	- else :
	- for ix in range(0,len(couplings_VVS)) :
	- if(tc[ix] == 0.) :
	- continue
	- elif(couplings_VVS[ix]==0.) :
	- couplings_VVS[ix]=tc[ix]
	- else :
	- couplings_VVS[ix] = '(( %s ) + ( %s ) )' % (couplings_VVS[ix],tc[ix])
	+
	+ elif lorentztag == 'VVS' or lorentztag == "VVSS" or lorentztag == "VSS" :
	+ all_couplings[color_idx] = scalarVectorCouplings(vertex,value,prefactors,L,lorentztag,pos,
	+ all_couplings[color_idx],append,kinematics)
	elif lorentztag == "SSS" or lorentztag == "SSSS" :
	prepend, kinematics, all_couplings[color_idx] = scalarCouplings(vertex,value,prefactors,L,lorentztag,pos,
	all_couplings[color_idx],prepend,kinematics)
	else:
	- if lorentztag == 'VSS':
	- if L.structure == 'P(1,3) - P(1,2)':
	- prefactors += ' * (-1)'
	- append = 'if(p2->id()!=%s){norm(-norm());}' \
	- % vertex.particles[1].pdg_code
	- elif lorentztag == 'VVVV':
	+ if lorentztag == 'VVVV':
	if qcd==2:
	append = 'setType(1);\nsetOrder(0,1,2,3);'
	else:
	append, factor = EWVVVVCouplings(vertex,L)
	prefactors += ' * (%s)' % factor
	elif lorentztag == 'VVV':
	if len(pos[8]) != 3:
	append = VVVordering(vertex)
	elif lorentztag == 'VVVS' :
	if len(pos[8]) == 0 :
	append = VVVordering(vertex)
	coup_norm.append('(%s) * (%s)' % (prefactors,value))

	- # final processing of the couplings
	- if('T' in lorentztag) :
	- (coup_left,coup_right,coup_norm) = processTensorCouplings(lorentztag,vertex,self.model,
	- self.parmsubs,all_couplings)
	- elif(lorentztag=="SSS" or lorentztag=="SSSS") :
	- coup_norm = processScalarCouplings(lorentztag,vertex,self.model,self.parmsubs,all_couplings)
	- # return the result
	- return (coup_left,coup_right,coup_norm,couplings_VVS,kinematics,qcd,qed,prepend,append)
	+ # special for vertices still in old structure till they get moved
	+ if( lorentztag=="FFV" or lorentztag == "FFS" or lorentztag=="VVV" or lorentztag =="VVVS" or lorentztag == "VVVV") :
	+ leftcontent = ' + '.join(coup_left) if coup_left else '0'
	+ rightcontent = ' + '.join(coup_right) if coup_right else '0'
	+ normcontent = ' + '.join(coup_norm) if coup_norm else '1'
	+ all_couplings=[]
	+ all_couplings.append(normcontent)
	+ all_couplings.append(leftcontent)
	+ all_couplings.append(rightcontent)
	+
	+ # return the result
	+ return (all_couplings,kinematics,qcd,qed,prepend,append)

File Metadata

Mime Type: text/x-diff
Expires: Thu, Apr 24, 6:32 AM (1 d, 12 h)
Storage Engine: blob
Storage Format: Raw Data
Storage Handle: 4872914
Default Alt Text: (59 KB)

No OneTemporaryActions

View Options

File Metadata

Event Timeline

No OneTemporary
Actions