diff --git a/.travis.yml b/.travis.yml deleted file mode 100644 index 286fd4a..0000000 --- a/.travis.yml +++ /dev/null @@ -1,27 +0,0 @@ -# Documentation: http://docs.travis-ci.com/user/languages/julia/ -language: julia -os: - - linux - - osx -julia: - - 1.0 - - 1.1 - - 1.2 - - 1.3 - - 1.4 - - 1.5 - - nightly -notifications: - email: false -after_success: - - julia -e 'using Pkg; Pkg.add("Coverage"); using Coverage; Codecov.submit(process_folder())' -jobs: - include: - - stage: Documentation - julia: 1.0 - script: julia --project=docs -e ' - using Pkg; - Pkg.develop(PackageSpec(path=pwd())); - Pkg.instantiate(); - include("docs/make.jl");' - after_success: skip diff --git a/docs/make.jl b/docs/make.jl index b3b31e7..c1d68e1 100644 --- a/docs/make.jl +++ b/docs/make.jl @@ -9,6 +9,7 @@ makedocs(; ), pages=[ "Introduction" => "index.md", + "Helper Functions" => "helpers.md", "Acknowledgements" => "acknowledgements.md", "API" => "api.md", ], diff --git a/docs/src/helpers.md b/docs/src/helpers.md new file mode 100644 index 0000000..8a291b7 --- /dev/null +++ b/docs/src/helpers.md @@ -0,0 +1,45 @@ +# Helper Functions + +There are many useful functions available to check additional properties of +particles, like their types, composition, charge etc. + +All of them take either a `Particle`, a `PDGID` or a simple `Integer` (or +anything which can be converted to an `Integer`) which represents a PDG ID, as +input, so that the general API is the following: + +```julia +helperfunction(p::Union{Particle, PDGID, Integer}) +``` + +Here is a list of the currently available helper functions: + +- `hasup(p)` +- `hasdown(p)` +- `hasstange(p)` +- `hascharm(p)` +- `hasbottom(p)` +- `hastop(p)` +- `isquark(p)` +- `isstandard(p)` +- `isfundamental(p)` +- `fundamentalid(p)` +- `islepton(p)` +- `ismeson(p)` +- `isbaryon(p)` +- `ishadron(p)` +- `isRhadron(p)` +- `isSUSY(p)` +- `ispentaquark(p)` +- `isgaugebosonorhiggs(p)` +- `issmgaugebosonorhiggs(p)` +- `istechnicolor(p)` +- `iscompositequarkorlepton(p)` +- `isdyon(p)` +- `isdiquark(p)` +- `isgeneratorspecific(p)` +- `isspecial(p)` +- `isQball(p)` +- `hasfundamentalanti(p)` +- `isnucleus(p)` +- `A(p)` +- `Z(p)` diff --git a/docs/src/index.md b/docs/src/index.md index 55e3183..8020ca1 100644 --- a/docs/src/index.md +++ b/docs/src/index.md @@ -67,6 +67,35 @@ julia> p.mass 2045.0 MeV ± 9.0 MeV ``` +There are tons of helper functions to check other properties: + +```julia +julia> filter(hasstrange, particles()) +257-element Array{Particle,1}: + Particle(3224) Sigma(1385) + Particle(23124) Lambda(1890) + Particle(-13324) Xi(1820) + Particle(-329) K(4)*(2045) + ⋮ + Particle(13124) Lambda(1690) + Particle(-100321) K(1460) + Particle(20433) D(s1)(2460) +julia> filter(islepton, particles()) +16-element Array{Particle,1}: + Particle(18) nu(tau') + Particle(15) tau + Particle(11) e + Particle(13) mu + Particle(14) nu(mu) + Particle(17) tau' + ⋮ + Particle(12) nu(e) + Particle(-14) nu(mu) + Particle(-16) nu(tau) + Particle(-12) nu(e) + Particle(-13) mu +``` + ## Units For some properties like `mass` and `width` we use the diff --git a/src/Corpuscles.jl b/src/Corpuscles.jl index e19c958..a09a6de 100644 --- a/src/Corpuscles.jl +++ b/src/Corpuscles.jl @@ -10,8 +10,21 @@ import Base export Particle, PDGID, PythiaID, Geant3ID, particles +# helpers.jl +export isfundamental, isstandard +export isquark, islepton, ismeson, isbaryon, ishadron +export isRhadron, isSUSY, ispentaquark, isdyon, isnucleus, isdiquark +export istechnicolor, iscompositequarkorlepton +export isgaugebosonorhiggs, issmgaugebosonorhiggs +export isgeneratorspecific, isspecial, isQball, hasfundamentalanti +export hasdown, hasup, hascharm, hasstrange, hasbottom, hastop +export A, Z, charge, threecharge, J, S, L, jspin, sspin, lspin + # Julia 1.0 compatibility eachrow_(x) = (x[i, :] for i in 1:size(x)[1]) +isnothing(::Any) = false +isnothing(::Nothing) = true + const _data_dir = abspath(joinpath(@__DIR__, "..", "data")) @@ -139,6 +152,10 @@ struct Particle latex::String end +pdgid(p::PDGID) = p +pdgid(p::Particle) = p.pdgid +pdgid(x) = PDGID(Integer(x)) + function read_conversion_csv(filepath::AbstractString) file_content = readdlm(filepath, ',', AbstractString, skipstart=2, comments=true) conversions = parse.(Int, file_content[:,1:3]) @@ -293,6 +310,10 @@ function Base.show(io::IO, m::MeasuredValue) return end +function Base.show(io::IO, p::PDGID) + Printf.@printf(io, "PDGID(%s)", p.value) +end + function Base.show(io::IO, p::Particle) Printf.@printf(io, "Particle(%s) %s", p.pdgid.value, p.name) end @@ -317,4 +338,7 @@ function Base.print(io::IO, p::Particle) end end + +include("helpers.jl") + end # module diff --git a/src/helpers.jl b/src/helpers.jl new file mode 100644 index 0000000..0a0174f --- /dev/null +++ b/src/helpers.jl @@ -0,0 +1,731 @@ +# Helper functions +# The implementation is based on the `particle` Python package from +# the SciKit-HEP group: https://github.com/scikit-hep/particle with the Git +# commit 6fe8211181d6a0b21933b2d1042eeb4e87d77b3c from 2020-11-27 + +""" + isvalid(p::PDGID) + +Returns true if the PDG ID is valid. +""" +function Base.isvalid(p::PDGID) + fundamentalid(p) != 0 && return true + ismeson(p) && return true + isbaryon(p) && return true + isgaugebosonorhiggs(p) && return true + ispentaquark(p) && return true + isSUSY(p) && return true + isRhadron(p) && return true + isdyon(p) && return true + isdiquark(p) && return true + ispentaquark(p) && return true + isgeneratorspecific(p) && return true + istechnicolor(p) && return true + iscompositequarkorlepton(p) && return true + !isstandard(p) && return (isQball(p) || isnucleus(p)) + false +end + +""" + isstandard(p::Union{Particle, PDGID, Integer}) + +Returns true if the PDG ID of the particle follows the standard numbering scheme. +""" +function isstandard(p) + p = pdgid(p) + p.N8 == 0 && p.N9 == 0 && p.N10 == 0 +end + +""" + isfundamental(p::Union{Particle, PDGID, Integer}) + +""" +function isfundamental(p) + p = pdgid(p) + !isstandard(p) && return false + p.Nq2 == 0 && p.Nq1 == 0 && return true + abs(p.value) <= 100 && return true + false +end + + +""" + fundamentalid(p::Union{Particle, PDGID, Integer}) + +""" +function fundamentalid(p) + p = pdgid(p) + !isstandard(p) && return 0 + abspdgid = abs(p.value) + p.Nq2 == 0 && p.Nq1 == 0 && return abspdgid % 10000 + abs(p.value) <= 100 && return abspgdid + 0 +end + + +""" + isquark(p::Union{Particle, PDGID, Integer}) + +""" +isquark(p) = 1 <= abs(pdgid(p).value) <= 8 + +""" + islepton(p::Union{Particle, PDGID, Integer}) + +""" +function islepton(p) + p = pdgid(p) + !isstandard(p) && return false + 11 <= fundamentalid(p) <= 18 && return true + false +end + +""" + ismeson(p::Union{Particle, PDGID, Integer}) + +""" +function ismeson(p) + p = pdgid(p) + !isstandard(p) && return false + + abspdgid = abs(p.value) + + abspdgid <= 100 && return false + 0 < fundamentalid(p) <= 100 && return false + # Special IDs - K(L)0, ???, K(S)0 + abspdgid ∈ [130, 210, 310] && return true + # Special IDs - B(L)0, B(sL)0, B(H)0, B(sH)0 + abspdgid ∈ [150, 350, 510, 530] && return true + # Special particles - reggeon, pomeron, odderon + p.value ∈ [110, 990, 9990] && return true + + if p.Nj > 0 && p.Nq3 > 0 && p.Nq2 > 0 && p.Nq1 == 0 + # check for illegal antiparticles + p.Nq3 == p.Nq2 && p.value < 0 && return false + return true + end + false +end + +""" + isbaryon(p::Union{Particle, PDGID, Integer}) + +""" +function isbaryon(p) + p = pdgid(p) + abspdgid = abs(p.value) + + abspdgid <= 100 && return false + + # Special case of proton and neutron: + # needs to be checked first since isstandard() is false for nuclei + abspdgid ∈ [1000000010, 1000010010] && return true + + !isstandard(p) && return false + 0 < fundamentalid(p) <= 100 && return false + + # Old codes for diffractive p and n (MC usage) + abspdgid ∈ [2110, 2210] && return true + + p.Nj > 0 && p.Nq3 > 0 && p.Nq2 > 0 && p.Nq1 > 0 && return true + + (isRhadron(p) || ispentaquark(p)) && return false + + false +end + +""" + ishadron(p::Union{Particle, PDGID, Integer}) + +""" +function ishadron(p) + p = pdgid(p) + abs(p.value) ∈ [1000000010, 1000010010] && return true + !isstandard(p) && return false + ismeson(p) && return true + isbaryon(p) && return true + ispentaquark(p) && return true + isRhadron(p) && return true + false +end + +""" + isRhadron(p::Union{Particle, PDGID, Integer}) + +An R-hadron is of the form 10abcdj, 100abcj, or 1000abj, +where j = 2J + 1 gives the spin; b, c, and d are quarks or gluons; +and a (the digit following the zero's) is a SUSY particle. + +""" +function isRhadron(p) + p = pdgid(p) + !isstandard(p) && return false + (p.N != 1 || p.Nr != 0) && return false + isSUSY(p) && return false + # All R-hadrons have at least 3 core digits + (p.Nq2 == 0 || p.Nq3 == 0 || p.Nj == 0) && return false + true +end + + +""" + isSUSY(p::Union{Particle, PDGID, Integer}) + +Fundamental SUSY particles have N = 1 or 2. + +""" +function isSUSY(p) + p = pdgid(p) + !isstandard(p) && return false + p.N != 1 && p.N != 2 && return false + p.Nr != 0 && return false + fundamentalid(p) == 0 && return false + true +end + + +""" + ispentaquark(p::Union{Particle, PDGID, Integer}) + +Pentaquark IDs are of the form +/- 9 Nr Nl Nq1 Nq2 Nq3 Nj, +where Nj = 2J + 1 gives the spin and Nr Nl Nq1 Nq2 Nq3 denote the quark numbers +in order Nr >= Nl >= Nq1 >= Nq2 and Nq3 gives the antiquark number. + +""" +function ispentaquark(p) + p = pdgid(p) + !isstandard(p) && return false + p.N != 9 && return false + (p.Nr == 9 || p.Nr == 0) && return false + (p.Nj == 9 || p.Nl == 0) && return false + p.Nq1 == 0 && return false + p.Nq2 == 0 && return false + p.Nq3 == 0 && return false + p.Nj == 0 && return false + p.Nq2 > p.Nq1 && return false + p.Nq1 > p.Nl && return false + p.Nl > p.Nr && return false + true +end + +""" + isgaugebosonorhiggs(p::Union{Particle, PDGID, Integer}) + +""" +function isgaugebosonorhiggs(p) + p = pdgid(p) + 21 <= abs(p.value) <= 40 +end + +""" + issmgaugebosonorhiggs(p::Union{Particle, PDGID, Integer}) +""" +function issmgaugebosonorhiggs(p) + p = pdgid(p) + abs(p.value) == 24 || 21 <= p.value <= 25 +end + + +""" + istechnicolor(p::Union{Particle, PDGID, Integer}) +""" +function istechnicolor(p) + p = pdgid(p) + !isstandard(p) && return false + p.N == 3 +end + +""" + iscompositequarkorlepton(p::Union{Particle, PDGID, Integer}) + +Excited (composite) quarks and leptons have N = 4 and Nr = 0. +""" +function iscompositequarkorlepton(p) + p = pdgid(p) + !isstandard(p) && return false + fundamentalid(p) == 0 && return false + !(p.N == 4 && p.Nr == 0) && return false + true +end + +""" + isdyon(p::Union{Particle, PDGID, Integer}) + +Magnetic monopoles and Dyons are assumed to have one unit of Dirac monopole +charge and a variable integer number xyz units of electric charge, where xyz +stands for Nq1 Nq2 Nq3. +Codes 411xyz0 are used when the magnetic and electrical charge sign agree and +412xyz0 when they disagree, with the overall sign of the particle set by the +magnetic charge. For now, no spin information is provided. +""" +function isdyon(p) + p = pdgid(p) + !isstandard(p) && return false + p.N != 4 && return false + p.Nr != 1 && return false + !(p.Nl ∈ [1, 2]) && return false + p.Nq3 == 0 && return false + p.Nj != 0 && return false + true +end + +""" + isdiquark(p::Union{Particle, PDGID, Integer}) +""" +function isdiquark(p) + p = pdgid(p) + !isstandard(p) && return false + abs(p.value) <= 100 && return false + 0 < fundamentalid(p) <= 100 && return false + p.Nj > 0 && p.Nq3 == 0 && p.Nq2 > 0 && p.Nq1 > 0 && return true + false +end + +""" + isgeneratorspecific(p::Union{Particle, PDGID, Integer}) + +Codes 81-100 are reserved for generator-specific pseudoparticles and concepts. +Codes 901-930, 1901-1930, 2901-2930, and 3901-3930 are for additional components +of Standard Model parton distribution functions, where the latter three ranges +are intended to distinguish left/right/longitudinal components. +Codes 998 and 999 are reserved for GEANT tracking purposes. +""" +function isgeneratorspecific(p) + p = pdgid(p) + abspdgid = abs(p.value) + 81 <= abspdgid <= 100 && return true + 901 <= abspdgid <= 930 && return true + 1901 <= abspdgid <= 1930 && return true + 2901 <= abspdgid <= 2930 && return true + 3901 <= abspdgid <= 3930 && return true + abspdgid ∈ [998, 999] && return true + abspdgid ∈ [20022, 480000000] && return true # Special cases of opticalphoton and geantino + false +end + + +""" + isspecial(p::Union{Particle, PDGID, Integer}) + +Special particle in the sense of the classification in the PDG MC particle +numbering scheme document, hence the graviton, the DM (S = 0, 1/2, 1) particles, +the reggeons (reggeon, pomeron and odderon), and all generator-specific +pseudo-particles and concepts, see `isgeneratorspecific`. +""" +function isspecial(p) + p = pdgid(p) + p.value ∈ [39, 41, 42, 51, 52, 53, 110, 990, 9990] || isgeneratorspecific(p) +end + + +""" + isQball(p::Union{Particle, PDGID, Integer}) + +Does this PDG ID correspond to a Q-ball or any exotic particle with electric +charge beyond the qqq scheme? +Ad-hoc numbering for such particles is +/- 100XXXY0, where XXX.Y is the charge. +""" +function isQball(p) + p = pdgid(p) + p.N8 != 1 && return false + p.N != 0 && return false + p.Nr != 0 && return false + (abs(p.value) ÷ 10) % 10000 == 0 && return false + p.Nj != 0 && return false + true +end + + +""" + hasfundamentalanti(p::Union{Particle, PDGID, Integer}) + +If this is a fundamental particle, does it have a valid antiparticle? + +Notes +----- +Based on the current list of defined particles/concepts +in the PDG Monte Carlo Particle Numbering Scheme document. + +""" +function hasfundamentalanti(p) + p = pdgid(p) + fid = fundamentalid(p) # always a positive integer + + # Check generator-specific PDGIDs + 81 <= fid <= 100 && return fid ∈ [82, 84, 85, 86, 87] + + # Check PDGIDs from 1 to 79 + cp_conjugates = [21, 22, 23, 25, 32, 33, 35, 36, 39, 40, 43] + unassigned = vcat([9, 10, 19, 20, 26], 26:31, 45:79) # not in conversion.csv + + if (1 <= fid <= 79) && !(fid ∈ cp_conjugates) + fid ∈ unassigned && return false + return true + end + false +end + + +""" + isnucleus(p::Union{Particle, PDGID, Integer}) + +Ion numbers are +/- 10LZZZAAAI. +AAA is A - total baryon number +ZZZ is Z - total charge +L is the total number of strange quarks. +I is the isomer number, with I=0 corresponding to the ground state. +""" +function isnucleus(p) + p = pdgid(p) + # A proton can be a Hydrogen nucleus + # A neutron can be considered as a nucleus when given the PDG ID 1000000010, + # hence consistency demands that isnucleus(neutron) is True + abs(p.value) ∈ [2112, 2212] && return true + + if p.N10 == 1 && p.N9 == 0 + # Charge should always be less than or equal to the baryon number + A_ = A(p) + Z_ = Z(p) + + (isnothing(A_) || isnothing(Z_)) && return false + + A_ >= abs(Z_) && return true + end + false +end + + +""" + A(p::Union{Particle, PDGID, Integer}) + +Returns the atomic number A if the PDG ID corresponds to a nucleus, else it +returns `nothing`. +""" +function A(p) + p = pdgid(p) + abspdgid = abs(p.value) + abspdgid ∈ [2112, 2212] && return 1 + (p.N10 != 1 || p.N9 != 0) && return nothing + (abspdgid ÷ 10) % 1000 +end + + +""" + Z(p::Union{Particle, PDGID, Integer}) + +Returns the charge Z if the PDG ID corresponds to a nucleus, else it returns +`nothing`. +""" +function Z(p) + p = pdgid(p) + abspdgid = abs(p.value) + abspdgid == 2212 && return sign(p.value) + abspdgid == 2112 && return 0 + (p.N10 != 1 || p.N9 != 0) && return nothing + ((abspdgid ÷ 10000) % 1000) * sign(p.value) +end + +""" + jspin(p::Union{Particle, PDGID, Integer}) + +Returns the total spin as 2J+1. +""" +function jspin(p) + p = pdgid(p) + !isvalid(p) && return nothing + fid = fundamentalid(p) + if fid > 0 + 0 < fid < 7 && return 2 # 4th generation quarks not dealt with! + fid == 9 && return 3 + 10 < fid < 17 && return 2 # 4th generation quarks not deal with! + 20 < fid < 25 && return 3 # 4th generation quarks not deal with! + return nothing + end + abspdgid = abs(p.value) + abspdgid ∈ [1000000010, 1000010010] && return 2 # neutrion, proton + !isstandard(p) && return nothing + p.value ∈ [130, 310] && return 1 + return abspdgid % 10 +end + + +""" + J(p::Union{Particle, PDGID, Integer}) + +Returns the total spin J. +""" +function J(p) + jspin_ = jspin(p) + isnothing(jspin_) && return nothing + (jspin_ - 1) / 2 +end + + +""" + S(p::Union{Particle, PDGID, Integer}) + +Returns the spin S. + +This is valid for mesons only. `nothing` is returned otherwise. +Mesons with PDGIDs of the kind 9XXXXXX (N=9) are not experimentally well-known +particles and `nothing` is returned too. +""" +function S(p) + p = pdgid(p) + !ismeson(p) && return nothing + !isvalid(p) && return nothing + + abspdgid = abs(p.value) + (abspdgid ÷ 1000000) % 10 == 9 && return nothing # no knowledge so far + + nl = ((abspdgid) ÷ 10000) % 10 + js = abspdgid % 10 + + !(js == 1 || js >= 3) && return 0 + + if nl == 0 + js == 1 && return 0 + return 1 + end + if nl == 1 + js == 1 && return 1 + return 0 + end + if nl ∈ [2, 3] + js >= 3 && return 1 + return 0 + end + 0 +end + +""" + sspin(p::Union{Particle, PDGID, Integer}) + +Returns the spin S as 2S+1. + +This is valid for mesons only. `nothing` is returned otherwise. +Mesons with PDGIDs of the kind 9XXXXXX (N=9) are not experimentally well-known +particles and `nothing` is returned too. +""" +function sspin(p) + s = S(p) + isnothing(s) && return nothing + 2s + 1 +end + + +""" + L(p::Union{Particle, PDGID, Integer}) + +Returns the orbital angular momentum L. + +This is valid for mesons only. `nothing` is returned otherwise. +Mesons with PDGIDs of the kind 9XXXXXX (N=9) are not experimentally well-known +particles and `nothing` is returned too. +""" +function L(p) + p = pdgid(p) + !ismeson(p) && return nothing + !isvalid(p) && return nothing + + abspdgid = abs(p.value) + (abspdgid ÷ 1000000) % 10 == 9 && return nothing # no knowledge so far + + nl = (abspdgid ÷ 10000) % 10 + js = abspdgid % 10 + + if nl == 0 + js == 1 && return 0 + js == 3 && return 0 + js == 5 && return 1 + js == 7 && return 2 + js == 9 && return 3 + elseif nl == 1 + js == 1 && return 1 + js == 3 && return 1 + js == 5 && return 2 + js == 7 && return 3 + js == 9 && return 4 + elseif nl == 2 + js == 3 && return 1 + js == 5 && return 2 + js == 7 && return 3 + js == 9 && return 4 + elseif nl == 3 + js == 3 && return 2 + js == 5 && return 3 + js == 7 && return 4 + js == 9 && return 5 + end + + 0 +end + +""" + lspin(p::Union{Particle, PDGID, Integer}) + +Returns the orbital angular momentum L as 2S+1. + +This is valid for mesons only. `nothing` is returned otherwise. +Mesons with PDGIDs of the kind 9XXXXXX (N=9) are not experimentally well-known +particles and `nothing` is returned too. +""" +function lspin(p) + l = L(p) + isnothing(l) && return nothing + 2l + 1 +end + + + +""" + threecharge(p::Union{Particle, PDGID, Integer}) + +Returns 3 times the EM charge. +""" +function threecharge(p) + p = pdgid(p) + !isvalid(p) && return nothing + abspdgid = abs(p.value) + ch100 = [-1, 2, -1, 2, -1, 2, -1, 2, 0, 0, -3, 0, -3, 0, -3, 0, -3, 0, 0, 0, + 0, 0, 0, 3, 0, 0, 0, 0, 0, 0, 0, 0, 0, 3, 0, 0, 3, 0, 0, 0, 0, -1, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 6, 3, 6, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0] + + fid = fundamentalid(p) + + charge = nothing + + if !isstandard(p) + if isnucleus(p) # ion + Z_ = Z(p) + isnothing(Z_) && return nothing + return 3Z_ + elseif isQball(p) + charge = 3 * ((abspdgid ÷ 10) % 10000) + else + # this should never be reached in the present numbering scheme + # since extra bits exist only for Q-balls and nuclei + return nothing + end + elseif isdyon(p) + charge = 3 * ((abspdgid ÷ 10) % 1000) + # this is half right + # the charge sign will be changed below if pid < 0 + if p.Nl == 2 + charge = -charge + end + elseif 0 < fid <= 100 # use table from the SciKit-HEP particle project + charge = ch100[fid] + if abspdgid ∈ [1000017, 1000018, 1000034, 1000052, 1000053, 1000054] + charge = 0 + end + if abspdgid == 5100061 && abspdgid == 5100062 + charge = 6 + end + elseif p.Nj == 0 # KL, KS, or undefined + return 0 + elseif p.Nq1 == 0 || (isRhadron(p) && p.Nq1 == 9) # mesons + if p.Nq2 == 3 || p.Nq2 == 5 + charge = ch100[p.Nq3] - ch100[p.Nq2] + else + charge = ch100[p.Nq2] - ch100[p.Nq3] + end + elseif p.Nq3 == 0 # diquarks + charge = ch100[p.Nq2] + ch100[p.Nq1] + elseif isbaryon(p) || (isRhadron(p) && p.Nl == 9) # baryons + charge = ch100[p.Nq3] + ch100[p.Nq2] + ch100[p.Nq1] + end + + if !isnothing(charge) && p.value < 0 + charge = -charge + end + charge +end + +""" + charge(p::Union{Particle, PDGID, Integer}) + +Returns the EM charge. +""" +function charge(p) + c = threecharge(p) + !isnothing(c) && return c / 3 + return nothing +end + + +""" + _hasquark(p:::Union{Particle, PDGID, Integer}, q::Integer) + +Note that q is always positive, so [1, 6] for Standard Model quarks +and [7, 8] for fourth-generation quarks. +""" +function _hasquark(p, q::Integer) + p = pdgid(p) + # Nuclei can also contain strange quarks, + # cf. the definition of a nucleus PDG ID in isnucleus. + # This check needs to be done first since isstandard is false for nuclei + if isnucleus(p) + q ∈ [1, 2] && return true # Nuclei by construction contain up and down quarks + if q == 3 && !(p.value ∈ [2112, 2212]) + p.N8 > 0 && return true + return false + end + end + + !isstandard(p) && return false + fundamentalid(p) > 0 && return false + isdyon(p) && return false + + if isRhadron(p) + _digits = digits(abs(p.value), pad=10) + iz = 7 + for loc ∈ range(6; step=-1, stop=2) + if _digits[loc] == 0 + iz = loc + elseif loc == iz - 1 + # ignore squark or gluino + else + _digits[loc] == q && return true + end + end + return false + end + + (p.Nq3 == q || p.Nq2 == q || p.Nq1 == q) && return true + + ispentaquark(p) && (p.Nl == q || p.Nr == q) && return true + false +end + +""" + hasdown(p::Union{Particle, PDGID, Integer}) +""" +hasdown(p) = _hasquark(pdgid(p), 1) + +""" + hasdown(p::Union{Particle, PDGID, Integer}) +""" +hasup(p) = _hasquark(pdgid(p), 2) + +""" + hasstrange(p::Union{Particle, PDGID, Integer}) +""" +hasstrange(p) = _hasquark(pdgid(p), 3) + +""" + hascharm(p::Union{Particle, PDGID, Integer}) +""" +hascharm(p) = _hasquark(pdgid(p), 4) + +""" + hasbottom(p::Union{Particle, PDGID, Integer}) +""" +hasbottom(p) = _hasquark(pdgid(p), 5) + +""" + hastop(p::Union{Particle, PDGID, Integer}) +""" +hastop(p) = _hasquark(pdgid(p), 6) diff --git a/test/particles.jl b/test/particles.jl new file mode 100644 index 0000000..fe8ff74 --- /dev/null +++ b/test/particles.jl @@ -0,0 +1,119 @@ +# Particle IDs used for testing. The definition is the same as in +# the Python package 'particle' from the SciKit-HEP group: +# https://github.com/scikit-hep/particle/blob/master/tests/conftest.py + +@enum PDGIDS begin + # Gauge and Higgs bosons + Gluon = 21 + Photon = 22 + Z0 = 23 + WMinus = -24 + HiggsBoson = 25 + ZPrime = 32 + # Charged leptons + Electron = 11 + Positron = -11 + Muon = 13 + AntiMuon = -13 + Tau = 15 + TauPrime = 17 + # Neutrinos + Nu_e = 12 + NuBar_tau = -16 + # Quarks + DQuark = 1 + UQuark = 2 + SQuark = 3 + CQuark = 4 + BQuark = 5 + TQuark = 6 + BPrimeQuark = 7 # 4th generation + TPrimeQuark = 8 + # Quarkonia + jpsi = 443 + psi_2S = 100443 + Upsilon_1S = 553 + Upsilon_4S = 300553 + # Light hadrons + Pi0 = 111 + PiPlus = 211 + eta = 221 + eta_prime = 331 + a_0_1450_plus = 10211 + KL = 130 + KS = 310 + KMinus = -321 + rho_770_minus = -213 + phi = 333 + omega = 223 + K1_1270_0 = 10313 + K1_1400_0 = 20313 + rho_1700_0 = 30113 + a2_1320_minus = -215 + omega_3_1670 = 227 + f_4_2300 = 9010229 # example of a not-well-known meson + Proton = 2212 + AntiNeutron = -2112 + Lambda = 3122 + Sigma0 = 3212 + SigmaPlus = 3222 + SigmaMinus = 3112 + Xi0 = 3322 + AntiXiMinus = -3312 + OmegaMinus = 3334 + # Charm hadrons + D0 = 421 + DPlus = 411 + DsPlus = 431 + LcPlus = 4122 + # Beauty hadrons + B0 = 511 + BPlus = 521 + Bs = 531 + BcPlus = 541 + Lb = 5122 + # Top hadrons + T0 = 621 + LtPlus = 6122 + # Special particles + Graviton = 39 + Reggeon = 110 + Pomeron = 990 + Odderon = 9990 + # Supersymmetric particles + Gluino = 1000021 + Gravitino = 1000039 + STildeL = 1000003 + CTildeR = 2000004 + # R-hadrons + RPlus_TTildeDbar = 1000612 + R0_GTildeG = 1000993 + RPlusPlus_GTildeUUU = 1092224 + # Q-balls + QBall1 = 10000150 + QBall2 = -10000200 + # Dyons + DyonSameMagElecChargeSign = 4110010 + DyonOppositeMagElecChargeSign = 4120010 + # Di-quarks + DD1 = 1103 + SD0 = 3101 + # Nuclei + HydrogenNucleus = 1000010010 + Carbon12 = 1000060120 + # Pentaquarks + AntiUCbarCUDPentaquark = -9422144 + # example of spin 3/2 u-cbar-c-u-d pentaquark decaying to J/psi proton + UCbarCUDPentaquark = 9422144 + # Technicolor + Pi0TC = 3000111 + PiMinusTC = -3000211 + # Composite quarks and leptons + UQuarkStar = 4000002 + AntiElectronStar = -4000011 + # Generator specific pseudoparticles or concepts + AntiCHadron = -84 + # Invalid ID + Invalid1 = 0 # illegal ID + Invalid2 = 99999999 # general form is a 7-digit number +end diff --git a/test/runtests.jl b/test/runtests.jl index d5edc23..a26ad00 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -1,3 +1,6 @@ +# The tests are based on the test suite of the `particle` Python package from +# the SciKit-HEP group: https://github.com/scikit-hep/particle with the Git +# commit 6fe8211181d6a0b21933b2d1042eeb4e87d77b3c from 2020-11-27 using Test using InteractiveUtils using Corpuscles @@ -5,6 +8,9 @@ using Unitful const DATA_DIR = joinpath(@__DIR__, "../data") +include("particles.jl") + + @testset "inits" begin for T in vcat(map(subtypes, [Signed, Unsigned])...) p = Particle(T(11)) @@ -145,3 +151,395 @@ end @test occursin(r"PDG ID:\s*1", output) @test occursin(r"Name:\s*d\n", output) end + + +@testset "helpers" begin + Corpuscles.use_catalog_file(joinpath(DATA_DIR, "particle2020.csv")) + + + """ + Runs the given is/has-function on a set of candidates and also checks every + other noncandidate in PDGIDs for the opposite outcome. + """ + function check_candidates(f, candidates) + noncandidates = setdiff(Set(instances(PDGIDS)), Set(candidates)) + for candidate ∈ candidates + @test f(candidate) + end + for noncandidate ∈ noncandidates + @test !f(noncandidate) + end + end + + + for _particles ∈ [particles(), [p.pdgid for p in particles()], [p.pdgid.value for p in particles()]] + @test sum(map(isstandard, _particles)) > 0 + @test sum(map(isfundamental, _particles)) > 0 + @test sum(map(Corpuscles.isstandard, _particles)) == 610 + @test sum(map(Corpuscles.fundamentalid, _particles)) == 413 + + @test sum(map(isquark, _particles)) == 12 + @test sum(map(islepton, _particles)) == 16 + @test sum(map(ismeson, _particles)) == 234 + @test sum(map(isbaryon, _particles)) == 292 + @test sum(map(ishadron, _particles)) == 526 + @test sum(map(isRhadron, _particles)) == 1 + @test sum(map(isSUSY, _particles)) == 0 + @test sum(map(ispentaquark, _particles)) == 0 + @test sum(map(isgaugebosonorhiggs, _particles)) == 6 + @test sum(map(issmgaugebosonorhiggs, _particles)) == 6 + @test sum(map(isdyon, _particles)) == 0 + @test sum(map(isnucleus, _particles)) == 4 + @test sum(map(isdiquark, _particles)) == 50 + @test sum(map(istechnicolor, _particles)) == 0 + @test sum(map(iscompositequarkorlepton, _particles)) == 0 + @test sum(map(isgeneratorspecific, _particles)) == 0 + @test sum(map(isspecial, _particles)) == 0 + @test sum(map(isQball, _particles)) == 0 + @test sum(map(hasfundamentalanti, _particles)) == 30 + + @test sum(map(hasdown, _particles)) == 328 + @test sum(map(hasup, _particles)) == 346 + @test sum(map(hascharm, _particles)) == 107 + @test sum(map(hasstrange, _particles)) == 257 + @test sum(map(hasbottom, _particles)) == 68 + @test sum(map(hastop, _particles)) == 0 + end + + @testset "isquark" begin + candidates = [DQuark, UQuark, SQuark, CQuark, BQuark, TQuark, BPrimeQuark, TPrimeQuark] + check_candidates(isquark, candidates) + end + + @testset "islepton" begin + candidates = [Electron, Positron, Muon, AntiMuon, Tau, TauPrime, Nu_e, NuBar_tau, AntiElectronStar] + check_candidates(islepton, candidates) + end + + @testset "ismeson" begin + candidates = [jpsi, psi_2S, Upsilon_1S, Upsilon_4S, Pi0, PiPlus, eta, + eta_prime, a_0_1450_plus, KL, KS, KMinus, phi, omega, + rho_770_minus, K1_1270_0, K1_1400_0, rho_1700_0, a2_1320_minus, + omega_3_1670, f_4_2300, D0, DPlus, DsPlus, B0, BPlus, Bs, + BcPlus, Pi0TC, PiMinusTC, T0, Reggeon, Pomeron, Odderon, + RPlus_TTildeDbar, R0_GTildeG] + check_candidates(ismeson, candidates) + end + + @testset "isbaryon" begin + candidates = [Proton, AntiNeutron, HydrogenNucleus, Lambda, Sigma0, + SigmaPlus, SigmaMinus, Xi0, AntiXiMinus, OmegaMinus, + LcPlus, Lb, LtPlus, RPlusPlus_GTildeUUU, + UCbarCUDPentaquark, AntiUCbarCUDPentaquark] + check_candidates(isbaryon, candidates) + end + + @testset "ishadron" begin + @test all(ishadron(p) == (ismeson(p) || isbaryon(p)) for p in instances(PDGIDS)) + end + + @testset "ispentaquark" begin + candidates = [UCbarCUDPentaquark, AntiUCbarCUDPentaquark] + check_candidates(ispentaquark, candidates) + end + + @testset "isgaugebosonorhiggs" begin + candidates = [Gluon, Photon, Z0, WMinus, HiggsBoson, ZPrime, Graviton] + check_candidates(isgaugebosonorhiggs, candidates) + end + + @testset "issmgaugebosonorhiggs" begin + candidates = [Gluon, Photon, Z0, WMinus, HiggsBoson] + check_candidates(issmgaugebosonorhiggs, candidates) + end + + @testset "isgeneratorspecific" begin + candidates = [AntiCHadron] + check_candidates(isgeneratorspecific, candidates) + end + + @testset "isspecial" begin + candidates = [Graviton, Reggeon, Pomeron, Odderon, AntiCHadron] + check_candidates(isspecial, candidates) + end + + @testset "isnucleus" begin + candidates = [Proton, AntiNeutron, HydrogenNucleus, Carbon12] + check_candidates(isnucleus, candidates) + end + + @testset "isdiquark" begin + candidates = [DD1, SD0] + check_candidates(isdiquark, candidates) + end + + @testset "isRhadron" begin + candidates = [RPlus_TTildeDbar, R0_GTildeG, RPlusPlus_GTildeUUU] + check_candidates(isRhadron, candidates) + end + + @testset "isQball" begin + candidates = [QBall1, QBall2] + check_candidates(isQball, candidates) + end + + @testset "isdyon" begin + candidates = [DyonSameMagElecChargeSign, DyonOppositeMagElecChargeSign] + check_candidates(isdyon, candidates) + end + + @testset "isSUSY" begin + candidates = [Gluino, Gravitino, STildeL, CTildeR] + check_candidates(isSUSY, candidates) + end + + @testset "istechnicolor" begin + candidates = [Pi0TC, PiMinusTC] + check_candidates(istechnicolor, candidates) + end + + @testset "iscompositequarkorlepton" begin + candidates = [UQuarkStar, AntiElectronStar] + check_candidates(iscompositequarkorlepton, candidates) + end + + @testset "hasdown" begin + f = hasdown + @test f(Photon) == false + @test f(Gluon) == false + @test f(Electron) == false + @test f(AntiMuon) == false + @test f(jpsi) == false + @test f(Upsilon_1S) == false + @test f(PiPlus) == true + @test f(KMinus) == false + @test f(D0) == false + @test f(DPlus) == true + @test f(DsPlus) == false + @test f(B0) == true + @test f(Bs) == false + @test f(BcPlus) == false + @test f(Proton) == true + @test f(LcPlus) == true + @test f(Lb) == true + @test f(DD1) == true + @test f(SD0) == true + @test f(Invalid1) == false + @test f(Invalid2) == false + end + + @testset "hasup" begin + f = hasup + @test f(Photon) == false + @test f(Gluon) == false + @test f(Electron) == false + @test f(AntiMuon) == false + @test f(jpsi) == false + @test f(Upsilon_1S) == false + @test f(PiPlus) == true + @test f(KMinus) == true + @test f(D0) == true + @test f(DPlus) == false + @test f(DsPlus) == false + @test f(B0) == false + @test f(Bs) == false + @test f(BcPlus) == false + @test f(Proton) == true + @test f(LcPlus) == true + @test f(Lb) == true + @test f(DD1) == false + @test f(SD0) == false + @test f(Invalid1) == false + @test f(Invalid2) == false + end + + @testset "hasstrange" begin + f = hasstrange + @test f(Photon) == false + @test f(Gluon) == false + @test f(Electron) == false + @test f(AntiMuon) == false + @test f(jpsi) == false + @test f(Upsilon_1S) == false + @test f(PiPlus) == false + @test f(KMinus) == true + @test f(D0) == false + @test f(DPlus) == false + @test f(DsPlus) == true + @test f(B0) == false + @test f(Bs) == true + @test f(BcPlus) == false + @test f(Proton) == false + @test f(LcPlus) == false + @test f(Lb) == false + @test f(DD1) == false + @test f(SD0) == true + @test f(Invalid1) == false + @test f(Invalid2) == false + end + + @testset "hascharm" begin + candidates = [jpsi, psi_2S, D0, DPlus, DsPlus, BcPlus, LcPlus, UCbarCUDPentaquark, AntiUCbarCUDPentaquark] + check_candidates(hascharm, candidates) + end + + @testset "hasbottom" begin + candidates = [Upsilon_1S, Upsilon_4S, B0, BPlus, Bs, BcPlus, Lb] + check_candidates(hasbottom, candidates) + end + + @testset "hastop" begin + candidates = [T0, LtPlus] + check_candidates(hastop, candidates) + end + + @testset "hasfundamentalanti" begin + candidates = [WMinus, Electron, Positron, Muon, AntiMuon, Tau, TauPrime, + Nu_e, NuBar_tau, DQuark, UQuark, SQuark, CQuark, BQuark, + TQuark, BPrimeQuark, TPrimeQuark, UQuarkStar, + AntiElectronStar, STildeL, CTildeR, AntiCHadron] + check_candidates(hasfundamentalanti, candidates) + end + + @testset "A" begin + candidates = Dict(Proton => 1, AntiNeutron => 1, HydrogenNucleus => 1, Carbon12 => 12) + for (candidate, value) ∈ candidates + @test A(candidate) == value + end + noncandidates = setdiff(Set(keys(candidates)), Set(instances(PDGIDS))) + for noncandidate in noncandidates + @test isnothiung(A(candidate)) + end + end + + @testset "Z" begin + candidates = Dict(Proton => 1, AntiNeutron => 0, HydrogenNucleus => 1, Carbon12 => 6) + for (candidate, value) ∈ candidates + @test Z(candidate) == value + end + noncandidates = setdiff(Set(keys(candidates)), Set(instances(PDGIDS))) + for noncandidate in noncandidates + @test isnothiung(Z(candidate)) + end + end + + @testset "isvalid" begin + f = isvalid + candidates = [Photon, Gluon, Electron, AntiMuon, jpsi, Upsilon_1S, + PiPlus, KMinus, D0, DPlus, DsPlus, B0, Bs, BcPlus, Proton, + LcPlus, Lb, DD1, SD0] + noncandidates = [Invalid1, Invalid2] + for candidate ∈ candidates + @test f(Corpuscles.pdgid(candidate)) + end + for noncandidate in noncandidates + @test !f(Corpuscles.pdgid(noncandidate)) + end + end + + @testset "charge" begin + @test charge(Photon) == 0 + @test charge(Gluon) == 0 + @test charge(Electron) == -1 + @test charge(AntiMuon) == +1 + @test charge(jpsi) == 0 + @test charge(Upsilon_1S) == 0 + @test charge(PiPlus) == +1 + @test charge(KMinus) == -1 + @test charge(D0) == 0 + @test charge(DPlus) == +1 + @test charge(DsPlus) == +1 + @test charge(B0) == 0 + @test charge(Bs) == 0 + @test charge(BcPlus) == +1 + @test charge(Proton) == +1 + @test charge(LcPlus) == +1 + @test charge(Lb) == 0 + @test charge(DD1) == -2 / 3 + @test charge(SD0) == -2 / 3 + @test Corpuscles.isnothing(charge(Invalid1)) + @test Corpuscles.isnothing(charge(Invalid2)) + end + + @testset "threecharge" begin + @test threecharge(Photon) == 0 + @test threecharge(Electron) == -3 + @test threecharge(jpsi) == 0 + @test threecharge(Upsilon_1S) == 0 + @test threecharge(KMinus) == -3 + @test threecharge(D0) == 0 + @test threecharge(Proton) == +3 + @test threecharge(LcPlus) == +3 + @test threecharge(Lb) == 0 + @test threecharge(DD1) == -2 + @test threecharge(SD0) == -2 + @test Corpuscles.isnothing(threecharge(Invalid1)) + @test Corpuscles.isnothing(threecharge(Invalid2)) + end + + JSLstatelist = Dict( + "000" => [Pi0, PiPlus, eta, eta_prime, KL, KS, KMinus, D0, DPlus, + DsPlus, B0, BPlus, Bs, BcPlus, T0], + "011" => [a_0_1450_plus], + "101" => [K1_1270_0], + "110" => [rho_770_minus], + "111" => [K1_1400_0], + "112" => [rho_1700_0], + "211" => [a2_1320_minus], + "312" => [omega_3_1670], + ) + + @testset "JSL" begin + @testset "badly known meons" begin + @test jspin(f_4_2300) == 9 + @test sspin(f_4_2300) === nothing + @test sspin(f_4_2300) === nothing + end + + for (states, particles) ∈ JSLstatelist + for p in particles + J_, S_, L_ = [parse(Int, s) for s ∈ states] + @test J(p) == J_ + @test S(p) == S_ + @test L(p) == L_ + end + end + end + + + @testset "J non-mesons" begin + Jlist = Dict( + 1 => [Gluon, Photon, Z0, jpsi, psi_2S, Upsilon_1S, Upsilon_4S, K1_1270_0,], + 1//2 => [Electron, Positron, Muon, AntiMuon, Tau, Nu_e, NuBar_tau, + DQuark, UQuark, SQuark, CQuark, BQuark, TQuark, Proton, + AntiNeutron, Lambda, Sigma0, SigmaPlus, SigmaMinus, Xi0, + AntiXiMinus, LcPlus, Lb, LtPlus, STildeL, CTildeR], + 3//2 => [OmegaMinus], + nothing => [Invalid1, Invalid2], + nothing => [TauPrime, BPrimeQuark, TPrimeQuark] + ) + for (J_, particles) ∈ Jlist + for p in particles + @test J(p) == J_ + end + end + end + + @testset "S non-mesons" begin + Slist = Dict(nothing => [Gluon, Photon, Z0]) + for (S_, particles) ∈ Slist + for p in particles + @test S(p) == S_ + end + end + end + + @testset "L non-mesons" begin + Llist = Dict(nothing => [Gluon, Photon, Z0]) + for (L_, particles) ∈ Llist + for p in particles + @test L(p) == L_ + end + end + end +end