Skip to content

Commit

Permalink
Use faster findmin
Browse files Browse the repository at this point in the history
The fast version of findmin can also be used with the N2Plain
algorithm, giving a good speedup at higher densities
  • Loading branch information
graeme-a-stewart committed Nov 8, 2023
1 parent 47182d9 commit eb2466c
Show file tree
Hide file tree
Showing 29 changed files with 19 additions and 15,078 deletions.
2 changes: 1 addition & 1 deletion src/JetReconstruction.jl
Original file line number Diff line number Diff line change
Expand Up @@ -34,7 +34,7 @@ include("HepMC3.jl")
include("PlainAlgo.jl")
export plain_jet_reconstruct

## Tiled algorithms
## N2Tiled algorithm
# Common pieces
include("TiledAlgoUtils.jl")

Expand Down
10 changes: 3 additions & 7 deletions src/PlainAlgo.jl
Original file line number Diff line number Diff line change
@@ -1,3 +1,5 @@
using LoopVectorization

Base.@propagate_inbounds function dist(i, j, rapidity_array, phi_array)
drapidity = rapidity_array[i] - rapidity_array[j]
dphi = abs(phi_array[i] - phi_array[j])
Expand Down Expand Up @@ -104,7 +106,6 @@ function plain_jet_reconstruct(particles::Vector{T}; p = -1, R = 1.0, recombine
end



function plain_jet_reconstruct(;objects_array::Vector{J}, kt2_array::Vector{F},
phi_array::Vector{F}, rapidity_array::Vector{F}, p = -1, R = 1.0, recombine = +, ptmin = 0.0) where {J, F<:AbstractFloat}
# Bounds
Expand Down Expand Up @@ -137,12 +138,7 @@ function plain_jet_reconstruct(;objects_array::Vector{J}, kt2_array::Vector{F},
iteration::Int = 1
while N != 0
# findmin
i = 1
dij_min = nndij[1]
@inbounds @simd for k in 2:N
cond = nndij[k] < dij_min
dij_min, i = ifelse(cond, (nndij[k], k), (dij_min, i))
end
dij_min, i = fast_findmin(nndij, N)

j::Int = nn[i]

Expand Down
16 changes: 1 addition & 15 deletions src/TiledAlgoLL.jl
Original file line number Diff line number Diff line change
Expand Up @@ -259,20 +259,6 @@ function find_tile_neighbours!(tile_union, jetA, jetB, oldB, tiling)
n_near_tiles
end

"""Find the lowest value in the array, returning the value and the index"""
find_lowest(dij, n) = begin
best = 1
@inbounds dij_min = dij[1]
@turbo for here in 2:n
newmin = dij[here] < dij_min
best = newmin ? here : best
dij_min = newmin ? dij[here] : dij_min
end
# @assert argmin(dij[1:n]) == best
dij_min, best
end



"""Return all inclusive jets of a ClusterSequence with pt > ptmin"""
function inclusive_jets(clusterseq::ClusterSequence, ptmin = 0.0)
Expand Down Expand Up @@ -375,7 +361,7 @@ function tiled_jet_reconstruct(particles::Vector{PseudoJet}; p = -1, R = 1.0, re
# compact NNs and dij arrays
ilast = N - (iteration - 1)
# Search for the lowest value of min_dij_ijet
dij_min, ibest = find_lowest(dij, ilast)
dij_min, ibest = fast_findmin(dij, ilast)
@inbounds jetA = NNs[ibest]
jetB = jetA.NN

Expand Down
14 changes: 14 additions & 0 deletions src/Utils.jl
Original file line number Diff line number Diff line change
Expand Up @@ -127,3 +127,17 @@ function final_jets(jets::Vector{LorentzVectorCyl}, ptmin::AbstractFloat)
end
final_jets
end

## Faster findmin function
"""Find the lowest value in the array, returning the value and the index"""
fast_findmin(dij, n) = begin
# findmin(@inbounds @view dij[1:n])
best = 1
@inbounds dij_min = dij[1]
@turbo for here in 2:n
newmin = dij[here] < dij_min
best = newmin ? here : best
dij_min = newmin ? dij[here] : dij_min
end
dij_min, best
end
67 changes: 0 additions & 67 deletions test/data/1-fj-result.dat

This file was deleted.

Loading

0 comments on commit eb2466c

Please sign in to comment.