Possible bugs in Task 7

[ingomueller-net]

I have several doubts about the current implementation of Task 7.

The benchmark description says to plot "jets [...] that are not within 0.4 in ΔR of any lepton [...]", i.e., of any muon or electron. In math, this is equivalent to the following conditions, which are all equivalent among themselves (omitting the conditions on pt):

(1) {j \in Jets | \nexists l \in Leptons :    DeltaR(j, l) < 0.4 }
(2) {j \in Jets | (\nexists m \in Muons :     DeltaR(j, m) < 0.4) \and
                  (\nexists e \in Electrons : DeltaR(j, e) < 0.4) }
(3) {j \in Jets | (\forall m \in Muons :    DeltaR(j, m) >= 0.4) \and
                  (\forall e \in Electrons: DeltaR(j, e) >= 0.4) }
(4) {j \in Jets | \forall m \in Muons :     DeltaR(j, m) >= 0.4 } \intersect
    {j \in Jets | \forall e \in Electrons:  DeltaR(j, e) >= 0.4 }

The current implementations is based on goodJets, which, roughly speaking, implements the following (it is called once on with muons, once with electrons):

(5) { j \in Jets | \exists m \in Muons :    DeltaR(j, m) >= 0.4) }
(6) { j \in Jets | \exists e \in Electrons: DeltaR(j, e) >= 0.4) }

That isn't the same. (5) and (6) compute if there exists a good lepton; (1) through (4) compute if there exists no bad one.

I see two other problems: First, goodJets is used the following way:

var (
	jets    []int
	muJets  = goodJets(jetPt, jetEta, jetPhi, muPt, muEta, muPhi)
	eleJets = goodJets(jetPt, jetEta, jetPhi, elePt, eleEta, elePhi)
)
switch {
case len(muJets) > 0:  // If there is one good jet as far as muons are concerned, electrons do not matter
	jets = muJets
case len(eleJets) > 0:
	jets = eleJets
default:
	continue
}

The problem is that the two types of leptons are never looked at jointly (like in formulation (3) above) and the results testing for each type are in no way combined (like in formulation (4) above). It seems like one of the two is necessary to compute the correct result.

Second, goodMuon may return the index of the same jet several times. This is the relevant part:

var jets = make([]int, 0, len(pt1))

for ijet, jetPt := range pt1 {
	// ...
	for ilep, lepPt := range pt2 {
		// ...
		dr2 := ...
		if dr2 > dr2Min {
			jets = append(jets, ijet)
			// After we added ijet, we continue with the next lepton,
			// which might meet the condition again,
			// in which case ijet is added again
		}
	}
},

[ingomueller-net]

I have implemented a version of which I think it computes the desired result. I donate it to the public domain should you want to use it. The main part of basic7 should look like this:

var (
    badMuJets  = badJets(jetPt, jetEta, jetPhi, muPt, muEta, muPhi)
    badEleJets = badJets(jetPt, jetEta, jetPhi, elePt, eleEta, elePhi)
)
badJets1 := make(map[int]struct{})
for _, ijet := range append(badMuJets, badEleJets...) {
    badJets1[ijet] = struct{}{}
}
fmt.Printf("nBadJets: %v, nJets: %v\n", len(badJets1), njets)

pt := 0.0
for i, _ := range jetPt {
    _, isBad := badJets1[i]
    // Skip over problematic jets
    if !isBad {
        pt += float64(jetPt[i])
    }
}
h1.Fill(pt, 1)

And we need this function:

func badJets(pt1, eta1, phi1, pt2, eta2, phi2 []float32) []int {
    const (
        twopi  = 2 * math.Pi
        dr2Min = 0.4 * 0.4
    )

    var jets = make([]int, 0, len(pt1))

    for ijet, jetPt := range pt1 {
        if jetPt <= 30 {
            // Not a jet with (pt > 30) --> we are not interested in it
            jets = append(jets, ijet)
            continue
        }
        // Search for bad leptons
        isBad := false
        for ilep, lepPt := range pt2 {
            if lepPt <= 10 {
                // Not a lepton with (pt > 10) --> not a problem
                continue
            }
            dphi := -math.Remainder(float64(phi1[ijet]-phi2[ilep]), twopi)
            deta := float64(eta1[ijet] - eta2[ilep])
            dr2 := dphi*dphi + deta*deta
            if dr2 < dr2Min {
                // Lepton within (DeltaR < 0.4) --> this is a bad jet
                isBad = true
            }
        }
        if isBad {
            // We found a bad lepton --> this jet is bad
            jets = append(jets, ijet)
        }
    }

    return jets
}

Let me know if you want me to create a PR.

[sbinet]

thanks for double-checking this :)

here is my take on it:

• groot/bench-opendata: fix task-7 #5