(* ::Package:: *)

(* MB integrand of massless box with one off-shell leg with p1^2 = M2 *)
OSBint = (-M2)^z2*(-t)^z4/(-s)^(z2+z4)*Gamma[2+eps+z2+z4]*Gamma[1+z2+z4]*Gamma[1+z4]*Gamma[-1-eps-z2-z4]*Gamma[-1-eps-z4]*Gamma[-z2]*Gamma[-z4];


(* Check that contours can be chosen straight for eps=-1, i.e. all Gamma functions have positive argument as ziAux\[Rule]0 *)
OSBintTilde = OSBint /. {eps->epsTil-1,z2->-1/2+z2Til,z4->-1/4+z4Til}


(* Collect terms as list *)
OSBintResult = {OSBintTilde};
(* Cross pole of Gamma[1/4-epsTil-z4Til] at z4Til = 1/4-epsTil *)
AppendTo[OSBintResult,-Residue[OSBintResult[[1]],{z4Til,1/4-epsTil}]]


(* Cross pole of Gamma[1/2-epsTil+z2Til] at z2Til = -1/2+epsTil *)
AppendTo[OSBintResult,+Residue[OSBintResult[[2]],{z2Til,-1/2+epsTil}]]


(* Cross pole of Gamma[3/4-epsTil-z2Til-z4Til] at z4Til = 3/4-epsTil-z2Til *)
AppendTo[OSBintResult,-Residue[OSBintResult[[1]],{z4Til,3/4-epsTil-z2Til}]]


(* Back to original variables with contours at Re[z2]=-1/2 and Re[z4]=-1/4 *)
OSBintResult = Simplify[OSBintResult/.{epsTil->1+eps,z2Til->1/2+z2,z4Til->1/4+z4}]


(* Pull out prefactor Exp[-eps*EulerGamma]/(st) in result *)
(* part without MB integral *)
Part1 = FullSimplify[Series[(s t Exp[eps*EulerGamma])/((-s)^(2+eps)Gamma[-2eps]) ((-M2)^eps s (-t)^-eps Gamma[-eps]^3 Gamma[1+eps]^2)/t,{eps,0,0}]]
(* part with 1 MB integral, close contour to the right *)
FullSimplify[-Residue[Normal[Series[OSBintResult[[2]]+OSBintResult[[4]],{eps,0,-1}]],{z2,n},Assumptions->n\[Element]Integers&&n>0]];
Part2 = FullSimplify[FullSimplify[Series[(s t Exp[eps*EulerGamma])/((-s)^(2+eps)Gamma[-2eps])*(-Residue[Normal[Series[OSBintResult[[2]]+OSBintResult[[4]],{eps,0,-1}]],{z2,0}]+Sum[%,{n,1,Infinity}]),{eps,0,0}]] /. {
   PolyLog[2,M2/s]->-PolyLog[2,1-M2/s]+\[Pi]^2/6-(Log[-M2]-Log[-s])Log[1-M2/s],PolyLog[2,M2/t]->-PolyLog[2,1-M2/t]+\[Pi]^2/6-(Log[-M2]-Log[-t])Log[1-M2/t]}]
(* part with 2 MB integrals is O[eps] *)


res = FullSimplify[Part1+Part2]


(* ::Subsubsection:: *)
(*Compare with Package-X*)


<<X`;
Series[s*t*
LoopRefine[LoopIntegrate[1,l,{l,0,1},{l+p1,0,1},{l+p1+p2,0,1},{l+p3,0,1}] /. {
p1.p1->M2,p2.p2->0,p3.p3->0,p1.p2->(s-M2)/2,p1.p3->(M2-t)/2,p2.p3->(s+t-M2)/2}] /. {
\[Micro]->1,\[Epsilon]->eps} /. {Log[-(1/M2)]->-Log[-M2],Log[-(1/s)]->-Log[-s],Log[-(1/t)]->-Log[-t],DiLog[x_,y_]->PolyLog[2,x]},
{eps,0,0}]
res - %


(* ::Subsubsection:: *)
(*Compare with MB.m*)


<<"~/applications/MB-1.2/MB.m"


OSB = OSBint/Gamma[-2*eps];


OSBrules = MBoptimizedRules[OSB,eps->0,{},{eps}]


OSBcont = MBcontinue[OSB,eps->0,OSBrules]


OSBselect = MBmerge[MBpreselect[OSBcont,{eps,0,0}]]


OSBexp = MBmerge[MBexpand[OSBselect,Exp[eps EulerGamma],{eps,0,0}]]