Subroutine ZCHARM ( kpi_pt , kpi_eta , & k3pi_pt , k3pi_eta, & D_s_pt , D_s_eta , & D0_pt , D0_eta , & kpi_pt_cr, kpi_eta_cr ) C# ---------------------------------------------------------------- C# C@# zcharm : ZES D* ( Kpi, K3pi ), D_s and D0 variables C# C# ---------------------------------------------------------------- C# C# input : ZEUS MDST C# C# ---------------------------------------------------------------- C# C# output : variables set to -1. in the case of no candidates C# C# kpi_pt, kpi_eta - Max(pt) and Min(|eta|) of all D* C# candidates with D0->Kpi C# ( CTD+REG ,right and wrong charges ) C# C# k3pi_pt, k3pi_eta - Max(pt) and Min(|eta|) of all D* C# candidates with D0->K3pi C# ( CTD+REG ,right and wrong charges ) C# C# D_s_pt, D_s_eta - Max(pt) and Min(|eta|) of all D_s C# candidates ( CTD+REG, right charges ) C# C# D0_pt, D0_eta - Max(pt) and Min(|eta|) of all D0 C# candidates ( CTD+REG, right charges ) C# C# kpi_pt_cr, kpi_eta_cr - Max(pt) and Min(|eta|) of all D* C# candidates with D0->Kpi C# ( CTD+REG , control region ) C# C# ---------------------------------------------------------------- C# Editor : Leonid Gladilin, gladilin@mail.desy.de C# C# Version 1.3, Date : 99-04-26 C# ---------------------------------------------------------------- * Implicit NONE * #include "partap.inc" #include "zdskey.inc" #include "zrevt.inc" #include "vctpar.inc" #include "vctvtx.inc" #include "vctrhl.inc" #include "vcevtctd.inc" #include "vcevtreg.inc" * Real kpi_pt, kpi_eta, k3pi_pt, k3pi_eta, D_s_pt, D_s_eta & , D0_pt , D0_eta , kpi_pt_cr, kpi_eta_cr * Integer n4v, i, j, k, j2, j3, jerr, ntrkus * Real piv , piby2, piml2 , doti, vdot , vmod, vdotn &, maspi / 0.13956995 /, dmaspi / 0.2791399 / &, sqmpi / 0.01947977 /, sqmka / 0.24371698 / &, gls(4), glt(4), vectus(4,300), rp1, rt1 &, vp1, vt1, vf1, vp2, vt2, vf2, vp3, vt3, vf3 &, vp4, vt4, vf4, vp5, vt5, vf5, vpt &, vr1(4), vr2(4), vr3(4), vr4(4), vr5(4) &, vmd0, ptd0, pmd0, etd0, tmd0 &, vmds, ptds, pmds, etds, tmds &, vmph, v2tmp, v3tmp, gtmp2(4), gtmp3(4) * Parameter ( piv = 3.1415926 ) Parameter ( piby2 = piv * 0.5 ) Parameter ( piml2 = piv * 2.0 ) * Integer ntrkm2 / 2 /, ntrkm3 / 3 /, ntrkm5 / 5 / * * Inlusive D0 selection criteria : * * only right charge combinations * * md0mni < M(D0) < md0mxi , Pt(tracks) > ptcuti, * Real md0mni / 1.59 /, md0mxi / 2.20 /, ptcuti / 0.45 / * * Inlusive D_s selection criteria : * * only right charge combinations * * ( M(Phi) < mphmx1 .and. mdsmn2 < M(D_s) < mdsmx2 ) .or. * ( mphmn2 < M(Phi) < mphmx2 .and. mdsmn1 < M(D_s) < mdsmx1 ) * & , mphmx1 / 1.06 / & , mphmn2 / 0.9 /, mphmx2 / 1.06 / & , mdsmn1 / 1.69 /, mdsmx1 / 2.2 / & , mdsmn2 / 1.78 /, mdsmx2 / 2.05 / * * Inlusive D* selection criteria : * * both right & wrong charge combinations ( and sum of charges = 1 ) * * ( delta_M < dmmax1 .and. md0mn1 < M(D0) < md0mx1 ) .or. * ( delta_M < dmmax2 .and. md0mn2 < M(D0) < md0mx2 ) * * control region for right charge combinations : * * ( delta_M < dmmax2 .and. md0mn3 < M(D0) < md0mx3 ) * & , md0mn1 / 1.39 /, md0mx1 / 2.20 / & , md0mn2 / 1.77 /, md0mx2 / 1.95 / & , md0mn3 / 1.95 /, md0mx3 / 2.20 / & , dmmax1 / 0.15 /, dmmax2 / 0.18 / * kpi_pt = -1. kpi_eta = -1. k3pi_pt = -1. k3pi_eta = -1. D_s_pt = -1. D_s_eta = -1. D0_pt = -1. D0_eta = -1. kpi_pt_cr = -1. kpi_eta_cr = -1. * if ( zdskey_gaftyp(1:3) .ne. 'EVT' ) Return * if ( coutab(VCEVTCTD) .le. 0 ) goto 999 * Call VCGetCTD(jerr) if ( jerr .ne. 0 ) then print 7,ZDSKEY_nr1,ZDSKEY_nr2,jerr 7 format(/,' !!! ZCHARM : ZDSKEY_nr1,nr2 =',i6,i7 &,' , VCGetCTD error =',i9,/) goto 999 endif * if ( COUTAB(VCTVTX) .le. 0 ) goto 999 * Call Fettab(VCTVTX,ID,1) * if ( VCTVTX_NDF .gt. 0 .and. VCTVTX_Chi2 .gt. 0. ) goto 9 print 8,ZDSKEY_nr1,ZDSKEY_nr2,VCTVTX_Chi2,VCTVTX_NDF 8 format(/,' !!! ZCHARM : corrupted CTD tracking : ZDSKEY_nr1 =' &,i6,' , ZDSKEY_nr2 =',i7, &/, ' VCTVTX_Chi2 =',e13.4,' , VCTVTX_NDF =' &,i4,/) goto 999 9 continue * n4v = COUTAB(VCTPAR) if ( n4v .lt. ntrkm2 ) goto 999 c ntrkus = 0 c do 11 i = 1 , n4v c Call Fettab(VCTPAR,ID,i) c Call Fettab(VCTRHL,ID,VCTPAR_VCTRHL) c rp1 = VCTRHL_pgevc rt1 = (piby2 - atan(vctrhl_tdip)) c vt1 = vctpar_par(1) vf1 = vctpar_par(2) if ( VCTPAR_par(3).eq.0. .or. sin(vt1).eq.0. ) goto 11 vpt = rp1*sin(rt1)*abs(VCTRHL_qovr/VCTPAR_par(3)) c vp1 = vpt / sin(vt1) c ntrkus = ntrkus + 1 vectus(1,ntrkus) = vt1 vectus(2,ntrkus) = vf1 vectus(3,ntrkus) = sign(vp1,VCTPAR_par(3)) vectus(4,ntrkus) = vpt c if ( ntrkus .eq. 300 ) goto 12 c 11 continue c if ( ntrkus .lt. ntrkm2 ) goto 999 c 12 do 301 i = 1 , ntrkus c vt1 = vectus(1,i) vf1 = vectus(2,i) vp1 = vectus(3,i) vpt = vectus(4,i) if ( vpt .lt. ptcuti ) goto 301 c vr1(1) = vpt* cos(vf1) vr1(2) = vpt* sin(vf1) vr1(3) = abs(vp1) * cos(vt1) vr1(4) = sqrt(vp1**2 + sqmka) c do 201 j = 1 , ntrkus if ( i .eq. j ) goto 201 c vt2 = vectus(1,j) vf2 = vectus(2,j) vp2 = vectus(3,j) vpt = vectus(4,j) if ( vpt .lt. ptcuti ) goto 201 c c remove wrong charge combinations c if ( vp1*vp2 .gt. 0. ) goto 201 c vr2(1) = vpt* cos(vf2) vr2(2) = vpt* sin(vf2) vr2(3) = abs(vp2) * cos(vt2) vr2(4) = sqrt(vp2**2 + sqmpi) c call vadd(vr1,vr2,gls,4) vmd0 = sqrt(-doti(gls,gls)) c if ( vmd0 .lt. md0mni .or. vmd0 .gt. md0mxi ) goto 201 c ptd0 = vmod(gls,2) pmd0 = vmod(gls,3) tmd0 = acos(gls(3)/pmd0) etd0 = abs(alog(tan(tmd0*0.5))) c c D0 candidates c if ( ptd0 .gt. D0_pt ) D0_pt = ptd0 if ( etd0 .lt. D0_eta .or. D0_eta .lt. 0. ) D0_eta = etd0 c 201 continue c 301 continue c if ( ntrkus .lt. ntrkm3 ) goto 999 c do 2301 i = 1 , ntrkus - 1 c vt1 = vectus(1,i) vf1 = vectus(2,i) vp1 = vectus(3,i) vpt = vectus(4,i) c vr1(1) = vpt* cos(vf1) vr1(2) = vpt* sin(vf1) vr1(3) = abs(vp1) * cos(vt1) vr1(4) = sqrt(vp1**2 + sqmka) c do 2201 j = i + 1 , ntrkus c vt2 = vectus(1,j) vf2 = vectus(2,j) vp2 = vectus(3,j) vpt = vectus(4,j) c c remove wrong charge combinations c if ( vp1*vp2 .gt. 0. ) goto 2201 c vr2(1) = vpt* cos(vf2) vr2(2) = vpt* sin(vf2) vr2(3) = abs(vp2) * cos(vt2) vr2(4) = sqrt(vp2**2 + sqmka) c call vadd(vr1,vr2,gls,4) vmph = sqrt(-doti(gls,gls)) c if ( vmph .gt. mphmx1 ) goto 2201 c do 2101 k = 1 , ntrkus if ( i .eq. k .or. & j .eq. k ) goto 2101 c vt3 = vectus(1,k) vf3 = vectus(2,k) vp3 = vectus(3,k) vpt = vectus(4,k) c vr3(1) = vpt* cos(vf3) vr3(2) = vpt* sin(vf3) vr3(3) = abs(vp3) * cos(vt3) vr3(4) = sqrt(vp3**2 + sqmpi) c call vadd(vr3,gls,glt,4) vmds = sqrt(-doti(glt,glt)) c if ( vmds .lt. mdsmn1 .or. vmds .gt. mdsmx1 ) goto 2101 c if ( ( vmph .lt. mphmn2 .or. vmph .gt. mphmx2 ) .and. & ( vmds .lt. mdsmn2 .or. vmds .gt. mdsmx2 ) ) goto 2101 c ptds = vmod(glt,2) pmds = vmod(glt,3) tmds = acos(glt(3)/pmds) etds = abs(alog(tan(tmds*0.5))) c c D_s candidates c if ( ptds .gt. D_s_pt ) D_s_pt = ptds if ( etds .lt. D_s_eta .or. D_s_eta .lt. 0. ) D_s_eta = etds c 2101 continue c 2201 continue c 2301 continue c do 3301 i = 1 , ntrkus c vt1 = vectus(1,i) vf1 = vectus(2,i) vp1 = vectus(3,i) vpt = vectus(4,i) c vr1(1) = vpt* cos(vf1) vr1(2) = vpt* sin(vf1) vr1(3) = abs(vp1) * cos(vt1) vr1(4) = sqrt(vp1**2 + sqmka) c do 3201 j = 1 , ntrkus if ( i .eq. j ) goto 3201 c vt2 = vectus(1,j) vf2 = vectus(2,j) vp2 = vectus(3,j) vpt = vectus(4,j) c vr2(1) = vpt* cos(vf2) vr2(2) = vpt* sin(vf2) vr2(3) = abs(vp2) * cos(vt2) vr2(4) = sqrt(vp2**2 + sqmpi) c call vadd(vr1,vr2,gls,4) vmd0 = sqrt(-doti(gls,gls)) c if ( ( vmd0.lt.md0mn1 .or. vmd0.gt.md0mx1 ) ) goto 3201 c do 3101 k = 1 , ntrkus if ( i .eq. k .or. & j .eq. k ) goto 3101 c vt3 = vectus(1,k) vf3 = vectus(2,k) vp3 = vectus(3,k) vpt = vectus(4,k) c c remove ((++)+), ((--)-), ((+-)+) and ((-+)-) combinations c if ( vp1*vp3 .gt. 0. ) goto 3101 c vr3(1) = vpt* cos(vf3) vr3(2) = vpt* sin(vf3) vr3(3) = abs(vp3) * cos(vt3) vr3(4) = sqrt(vp3**2 + sqmpi) c call vadd(vr3,gls,glt,4) c vmds = sqrt(-doti(glt,glt)) c if ( vmds-vmd0 .gt. dmmax2 ) goto 3101 c ptds = vmod(glt,2) pmds = vmod(glt,3) tmds = acos(glt(3)/pmds) etds = abs(alog(tan(tmds*0.5))) c if ( vmds-vmd0 .le. dmmax1 .or. & ( vmd0 .ge. md0mn2 .and. vmd0 .le. md0mx2 ) ) then c c right & wrong charge D* candidates c if ( ptds .gt. kpi_pt ) kpi_pt = ptds if ( etds .lt. kpi_eta .or. kpi_eta .lt. 0. ) kpi_eta = etds c endif c if ( vp1*vp2 .lt. 0. ) then if ( vmd0 .gt. md0mn3 .and. vmd0 .le. md0mx3 ) then c c control region right charge D* candidates c if ( ptds.gt.kpi_pt_cr ) kpi_pt_cr = ptds if ( etds.lt.kpi_eta_cr .or. kpi_eta_cr.lt.0. ) kpi_eta_cr = etds c endif endif c 3101 continue c 3201 continue c 3301 continue c if ( ntrkus .lt. ntrkm5 ) goto 999 c do 4301 i = 1 , ntrkus c vt1 = vectus(1,i) vf1 = vectus(2,i) vp1 = vectus(3,i) vpt = vectus(4,i) c vr1(1) = vpt* cos(vf1) vr1(2) = vpt* sin(vf1) vr1(3) = abs(vp1) * cos(vt1) vr1(4) = sqrt(vp1**2 + sqmka) c do 4201 j = 1 , ntrkus - 2 if ( i .eq. j ) goto 4201 c vt2 = vectus(1,j) vf2 = vectus(2,j) vp2 = vectus(3,j) vpt = vectus(4,j) c vr2(1) = vpt* cos(vf2) vr2(2) = vpt* sin(vf2) vr2(3) = abs(vp2) * cos(vt2) vr2(4) = sqrt(vp2**2 + sqmpi) c call vadd(vr1,vr2,gtmp2,4) v2tmp = sqrt(-doti(gtmp2,gtmp2)) c if ( v2tmp .gt. md0mx1 - dmaspi ) goto 4201 c do 4191 j2 = j + 1 , ntrkus - 1 if ( i .eq. j2 ) goto 4191 c vt3 = vectus(1,j2) vf3 = vectus(2,j2) vp3 = vectus(3,j2) vpt = vectus(4,j2) c vr3(1) = vpt* cos(vf3) vr3(2) = vpt* sin(vf3) vr3(3) = abs(vp3) * cos(vt3) vr3(4) = sqrt(vp3**2 + sqmpi) c call vadd(vr3,gtmp2,gtmp3,4) v3tmp = sqrt(-doti(gtmp3,gtmp3)) c if ( v3tmp .gt. md0mx1 - maspi ) goto 4191 c do 4181 j3 = j2 + 1 , ntrkus if ( i .eq. j3 ) goto 4181 c vt4 = vectus(1,j3) vf4 = vectus(2,j3) vp4 = vectus(3,j3) vpt = vectus(4,j3) c c remove (++++) and (----) combinations c if ( vp1*vp2.gt.0. .and. & vp1*vp3.gt.0. .and. & vp1*vp4.gt.0. ) goto 4181 c vr4(1) = vpt* cos(vf4) vr4(2) = vpt* sin(vf4) vr4(3) = abs(vp4) * cos(vt4) vr4(4) = sqrt(vp4**2 + sqmpi) c call vadd(vr4,gtmp3,gls,4) vmd0 = sqrt(-doti(gls,gls)) c if ( vmd0 .lt. md0mn1 .or. vmd0 .gt. md0mx1 ) goto 4181 c do 4101 k = 1 , ntrkus if ( i .eq. k .or. j .eq. k .or. & j2 .eq. k .or. j3 .eq. k ) goto 4101 c vt5 = vectus(1,k) vf5 = vectus(2,k) vp5 = vectus(3,k) vpt = vectus(4,k) c c remove ((+---)-) and ((-+++)+) combinations c if ( vp1*vp2.lt.0. .and. & vp1*vp3.lt.0. .and. & vp1*vp4.lt.0. .and. & vp1*vp5.lt.0. ) goto 4101 c c remove ((+++-)+) and ((---+)-) combinations c if ( vp1*vp2.gt.0. .and. & vp1*vp3.gt.0. .and. & vp1*vp4.lt.0. .and. & vp1*vp5.gt.0. ) goto 4101 c c remove ((++-+)+) and ((--+-)-) combinations c if ( vp1*vp2.gt.0. .and. & vp1*vp3.lt.0. .and. & vp1*vp4.gt.0. .and. & vp1*vp5.gt.0. ) goto 4101 c c remove ((+-++)+) and ((-+--)-) combinations c if ( vp1*vp2.lt.0. .and. & vp1*vp3.gt.0. .and. & vp1*vp4.gt.0. .and. & vp1*vp5.gt.0. ) goto 4101 c vr5(1) = vpt* cos(vf5) vr5(2) = vpt* sin(vf5) vr5(3) = abs(vp5) * cos(vt5) vr5(4) = sqrt(vp5**2 + sqmpi) c call vadd(vr5,gls,glt,4) vmds = sqrt(-doti(glt,glt)) c if ( vmds-vmd0 .gt. dmmax2 ) goto 4101 c ptds = vmod(glt,2) pmds = vmod(glt,3) tmds = acos(glt(3)/pmds) etds = abs(alog(tan(tmds*0.5))) c if ( vmds-vmd0 .le. dmmax1 .or. & ( vmd0 .ge. md0mn2 .and. vmd0 .le. md0mx2 ) ) then c c right & wrong charge D* candidates c if ( ptds .gt. k3pi_pt ) k3pi_pt = ptds if ( etds .lt. k3pi_eta .or. k3pi_eta .lt. 0. ) k3pi_eta = etds c endif c 4101 continue c 4181 continue 4191 continue 4201 continue c 4301 continue * 999 if ( coutab(VCEVTREG) .le. 0 ) Return * Call VCGetREG(jerr) if ( jerr .ne. 0 ) then print 17,ZDSKEY_nr1,ZDSKEY_nr2,jerr 17 format(/,' !!! ZCHARM : ZDSKEY_nr1,nr2 =',i6,i7 &,' , VCGetREG error =',i9,/) Return endif * if ( COUTAB(VCTVTX) .le. 0 ) Return * Call Fettab(VCTVTX,ID,1) * if ( VCTVTX_NDF .gt. 0 .and. VCTVTX_Chi2 .gt. 0. ) goto 19 print 18,ZDSKEY_nr1,ZDSKEY_nr2,VCTVTX_Chi2,VCTVTX_NDF 18 format(/,' !!! ZCHARM : corrupted REG tracking : ZDSKEY_nr1 =' &,i6,' , ZDSKEY_nr2 =',i7, &/, ' VCTVTX_Chi2 =',e13.4,' , VCTVTX_NDF =' &,i4,/) Return 19 continue * n4v = COUTAB(VCTPAR) if ( n4v .lt. ntrkm2 ) Return c ntrkus = 0 c do 21 i = 1 , n4v c Call Fettab(VCTPAR,ID,i) c Call Fettab(VCTRHL,ID,VCTPAR_VCTRHL) c rp1 = VCTRHL_pgevc rt1 = (piby2 - atan(vctrhl_tdip)) c vt1 = vctpar_par(1) vf1 = vctpar_par(2) if ( VCTPAR_par(3).eq.0. .or. sin(vt1).eq.0. ) goto 21 vpt = rp1*sin(rt1)*abs(VCTRHL_qovr/VCTPAR_par(3)) c vp1 = vpt / sin(vt1) c ntrkus = ntrkus + 1 vectus(1,ntrkus) = vt1 vectus(2,ntrkus) = vf1 vectus(3,ntrkus) = sign(vp1,VCTPAR_par(3)) vectus(4,ntrkus) = vpt c if ( ntrkus .eq. 300 ) goto 22 c 21 continue c if ( ntrkus .lt. ntrkm2 ) Return c 22 do 10301 i = 1 , ntrkus c vt1 = vectus(1,i) vf1 = vectus(2,i) vp1 = vectus(3,i) vpt = vectus(4,i) if ( vpt .lt. ptcuti ) goto 10301 c vr1(1) = vpt* cos(vf1) vr1(2) = vpt* sin(vf1) vr1(3) = abs(vp1) * cos(vt1) vr1(4) = sqrt(vp1**2 + sqmka) c do 10201 j = 1 , ntrkus if ( i .eq. j ) goto 10201 c vt2 = vectus(1,j) vf2 = vectus(2,j) vp2 = vectus(3,j) vpt = vectus(4,j) if ( vpt .lt. ptcuti ) goto 10201 c c remove wrong charge combinations c if ( vp1*vp2 .gt. 0. ) goto 10201 c vr2(1) = vpt* cos(vf2) vr2(2) = vpt* sin(vf2) vr2(3) = abs(vp2) * cos(vt2) vr2(4) = sqrt(vp2**2 + sqmpi) c call vadd(vr1,vr2,gls,4) vmd0 = sqrt(-doti(gls,gls)) c if ( vmd0 .lt. md0mni .or. vmd0 .gt. md0mxi ) goto 10201 c ptd0 = vmod(gls,2) pmd0 = vmod(gls,3) tmd0 = acos(gls(3)/pmd0) etd0 = abs(alog(tan(tmd0*0.5))) c c D0 candidates c if ( ptd0 .gt. D0_pt ) D0_pt = ptd0 if ( etd0 .lt. D0_eta .or. D0_eta .lt. 0. ) D0_eta = etd0 c 10201 continue c 10301 continue c if ( ntrkus .lt. ntrkm3 ) Return c do 12301 i = 1 , ntrkus - 1 c vt1 = vectus(1,i) vf1 = vectus(2,i) vp1 = vectus(3,i) vpt = vectus(4,i) c vr1(1) = vpt* cos(vf1) vr1(2) = vpt* sin(vf1) vr1(3) = abs(vp1) * cos(vt1) vr1(4) = sqrt(vp1**2 + sqmka) c do 12201 j = i + 1 , ntrkus c vt2 = vectus(1,j) vf2 = vectus(2,j) vp2 = vectus(3,j) vpt = vectus(4,j) c c remove wrong charge combinations c if ( vp1*vp2 .gt. 0. ) goto 12201 c vr2(1) = vpt* cos(vf2) vr2(2) = vpt* sin(vf2) vr2(3) = abs(vp2) * cos(vt2) vr2(4) = sqrt(vp2**2 + sqmka) c call vadd(vr1,vr2,gls,4) vmph = sqrt(-doti(gls,gls)) c if ( vmph .gt. mphmx1 ) goto 12201 c do 12101 k = 1 , ntrkus if ( i .eq. k .or. & j .eq. k ) goto 12101 c vt3 = vectus(1,k) vf3 = vectus(2,k) vp3 = vectus(3,k) vpt = vectus(4,k) c vr3(1) = vpt* cos(vf3) vr3(2) = vpt* sin(vf3) vr3(3) = abs(vp3) * cos(vt3) vr3(4) = sqrt(vp3**2 + sqmpi) c call vadd(vr3,gls,glt,4) vmds = sqrt(-doti(glt,glt)) c if ( vmds .lt. mdsmn1 .or. vmds .gt. mdsmx1 ) goto 12101 c if ( ( vmph .lt. mphmn2 .or. vmph .gt. mphmx2 ) .and. & ( vmds .lt. mdsmn2 .or. vmds .gt. mdsmx2 ) ) goto 12101 c ptds = vmod(glt,2) pmds = vmod(glt,3) tmds = acos(glt(3)/pmds) etds = abs(alog(tan(tmds*0.5))) c c D_s candidates c if ( ptds .gt. D_s_pt ) D_s_pt = ptds if ( etds .lt. D_s_eta .or. D_s_eta .lt. 0. ) D_s_eta = etds c 12101 continue c 12201 continue c 12301 continue c do 13301 i = 1 , ntrkus c vt1 = vectus(1,i) vf1 = vectus(2,i) vp1 = vectus(3,i) vpt = vectus(4,i) c vr1(1) = vpt* cos(vf1) vr1(2) = vpt* sin(vf1) vr1(3) = abs(vp1) * cos(vt1) vr1(4) = sqrt(vp1**2 + sqmka) c do 13201 j = 1 , ntrkus if ( i .eq. j ) goto 13201 c vt2 = vectus(1,j) vf2 = vectus(2,j) vp2 = vectus(3,j) vpt = vectus(4,j) c vr2(1) = vpt* cos(vf2) vr2(2) = vpt* sin(vf2) vr2(3) = abs(vp2) * cos(vt2) vr2(4) = sqrt(vp2**2 + sqmpi) c call vadd(vr1,vr2,gls,4) vmd0 = sqrt(-doti(gls,gls)) c if ( ( vmd0.lt.md0mn1 .or. vmd0.gt.md0mx1 ) ) goto 13201 c do 13101 k = 1 , ntrkus if ( i .eq. k .or. & j .eq. k ) goto 13101 c vt3 = vectus(1,k) vf3 = vectus(2,k) vp3 = vectus(3,k) vpt = vectus(4,k) c c remove ((++)+), ((--)-), ((+-)+) and ((-+)-) combinations c if ( vp1*vp3 .gt. 0. ) goto 13101 c vr3(1) = vpt* cos(vf3) vr3(2) = vpt* sin(vf3) vr3(3) = abs(vp3) * cos(vt3) vr3(4) = sqrt(vp3**2 + sqmpi) c call vadd(vr3,gls,glt,4) c vmds = sqrt(-doti(glt,glt)) c if ( vmds-vmd0 .gt. dmmax2 ) goto 13101 c ptds = vmod(glt,2) pmds = vmod(glt,3) tmds = acos(glt(3)/pmds) etds = abs(alog(tan(tmds*0.5))) c if ( vmds-vmd0 .le. dmmax1 .or. & ( vmd0 .ge. md0mn2 .and. vmd0 .le. md0mx2 ) ) then c c right & wrong charge D* candidates c if ( ptds .gt. kpi_pt ) kpi_pt = ptds if ( etds .lt. kpi_eta .or. kpi_eta .lt. 0. ) kpi_eta = etds c endif c if ( vp1*vp2 .lt. 0. ) then if ( vmd0 .gt. md0mn3 .and. vmd0 .le. md0mx3 ) then c c control region for right charge D* candidates c if ( ptds.gt.kpi_pt_cr ) kpi_pt_cr = ptds if ( etds.lt.kpi_eta_cr .or. kpi_eta_cr.lt.0. ) kpi_eta_cr = etds c endif endif c 13101 continue c 13201 continue c 13301 continue c if ( ntrkus .lt. ntrkm5 ) Return c do 14301 i = 1 , ntrkus c vt1 = vectus(1,i) vf1 = vectus(2,i) vp1 = vectus(3,i) vpt = vectus(4,i) c vr1(1) = vpt* cos(vf1) vr1(2) = vpt* sin(vf1) vr1(3) = abs(vp1) * cos(vt1) vr1(4) = sqrt(vp1**2 + sqmka) c do 14201 j = 1 , ntrkus - 2 if ( i .eq. j ) goto 14201 c vt2 = vectus(1,j) vf2 = vectus(2,j) vp2 = vectus(3,j) vpt = vectus(4,j) c vr2(1) = vpt* cos(vf2) vr2(2) = vpt* sin(vf2) vr2(3) = abs(vp2) * cos(vt2) vr2(4) = sqrt(vp2**2 + sqmpi) c call vadd(vr1,vr2,gtmp2,4) v2tmp = sqrt(-doti(gtmp2,gtmp2)) c if ( v2tmp .gt. md0mx1 - dmaspi ) goto 14201 c do 14191 j2 = j + 1 , ntrkus - 1 if ( i .eq. j2 ) goto 14191 c vt3 = vectus(1,j2) vf3 = vectus(2,j2) vp3 = vectus(3,j2) vpt = vectus(4,j2) c vr3(1) = vpt* cos(vf3) vr3(2) = vpt* sin(vf3) vr3(3) = abs(vp3) * cos(vt3) vr3(4) = sqrt(vp3**2 + sqmpi) c call vadd(vr3,gtmp2,gtmp3,4) v3tmp = sqrt(-doti(gtmp3,gtmp3)) c if ( v3tmp .gt. md0mx1 - maspi ) goto 14191 c do 14181 j3 = j2 + 1 , ntrkus if ( i .eq. j3 ) goto 14181 c vt4 = vectus(1,j3) vf4 = vectus(2,j3) vp4 = vectus(3,j3) vpt = vectus(4,j3) c c remove (++++) and (----) combinations c if ( vp1*vp2.gt.0. .and. & vp1*vp3.gt.0. .and. & vp1*vp4.gt.0. ) goto 14181 c vr4(1) = vpt* cos(vf4) vr4(2) = vpt* sin(vf4) vr4(3) = abs(vp4) * cos(vt4) vr4(4) = sqrt(vp4**2 + sqmpi) c call vadd(vr4,gtmp3,gls,4) vmd0 = sqrt(-doti(gls,gls)) c if ( vmd0 .lt. md0mn1 .or. vmd0 .gt. md0mx1 ) goto 14181 c do 14101 k = 1 , ntrkus if ( i .eq. k .or. j .eq. k .or. & j2 .eq. k .or. j3 .eq. k ) goto 14101 c vt5 = vectus(1,k) vf5 = vectus(2,k) vp5 = vectus(3,k) vpt = vectus(4,k) c c remove ((+---)-) and ((-+++)+) combinations c if ( vp1*vp2.lt.0. .and. & vp1*vp3.lt.0. .and. & vp1*vp4.lt.0. .and. & vp1*vp5.lt.0. ) goto 14101 c c remove ((+++-)+) and ((---+)-) combinations c if ( vp1*vp2.gt.0. .and. & vp1*vp3.gt.0. .and. & vp1*vp4.lt.0. .and. & vp1*vp5.gt.0. ) goto 14101 c c remove ((++-+)+) and ((--+-)-) combinations c if ( vp1*vp2.gt.0. .and. & vp1*vp3.lt.0. .and. & vp1*vp4.gt.0. .and. & vp1*vp5.gt.0. ) goto 14101 c c remove ((+-++)+) and ((-+--)-) combinations c if ( vp1*vp2.lt.0. .and. & vp1*vp3.gt.0. .and. & vp1*vp4.gt.0. .and. & vp1*vp5.gt.0. ) goto 14101 c vr5(1) = vpt* cos(vf5) vr5(2) = vpt* sin(vf5) vr5(3) = abs(vp5) * cos(vt5) vr5(4) = sqrt(vp5**2 + sqmpi) c call vadd(vr5,gls,glt,4) vmds = sqrt(-doti(glt,glt)) c if ( vmds-vmd0 .gt. dmmax2 ) goto 14101 c ptds = vmod(glt,2) pmds = vmod(glt,3) tmds = acos(glt(3)/pmds) etds = abs(alog(tan(tmds*0.5))) c if ( vmds-vmd0 .le. dmmax1 .or. & ( vmd0 .ge. md0mn2 .and. vmd0 .le. md0mx2 ) ) then c c right & wrong charge D* candidates c if ( ptds .gt. k3pi_pt ) k3pi_pt = ptds if ( etds .lt. k3pi_eta .or. k3pi_eta .lt. 0. ) k3pi_eta = etds c endif c 14101 continue c 14181 continue 14191 continue 14201 continue c 14301 continue * End