Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

add a kerneldensity method #584

Open
wants to merge 5 commits into
base: main
Choose a base branch
from
Open

add a kerneldensity method #584

wants to merge 5 commits into from

Conversation

rafaqz
Copy link
Owner

@rafaqz rafaqz commented Dec 22, 2023

@anjelinejeline if you want to try this out you can check out this branch with ] add Rasters#kerneldensity and it should just work. You can try the code in the tests here to get started.

@anjelinejeline
Copy link

@rafaqz I have installed the extension as you suggested ] add Rasters#kerneldensity from the REPL

(@v1.9) pkg> add Rasters#kerneldensity
     Cloning git-repo `https://github.com/rafaqz/Rasters.jl.git`
    Updating git-repo `https://github.com/rafaqz/Rasters.jl.git`
   Resolving package versions...
    Updating `~/.julia/environments/v1.9/Project.toml`
  [a3a2b9e3] ~ Rasters v0.10.0 ⇒ v0.10.0 `https://github.com/rafaqz/Rasters.jl.git#kerneldensity`
    Updating `~/.julia/environments/v1.9/Manifest.toml`
  [a3a2b9e3] ~ Rasters v0.10.0 ⇒ v0.10.0 `https://github.com/rafaqz/Rasters.jl.git#kerneldensity`
Precompiling project...
  3 dependencies successfully precompiled in 15 seconds. 275 already precompiled.

(@v1.9) pkg>

but the function is not recognized

using Rasters

# Load the raster template 
ras=Rasters.Raster("raster.tif")

# Load the CSV file
df = CSV.read("df.csv", DataFrame)

# Discard the first column
df = select(df, Not(1))

# Extract x and y columns
x_values = df.x
y_values = df.y

# Extract x and y coordinates
coordinates = hcat(x_values, y_values)

# Compute the kernel 

try1=Rasters.kerneldensity(data=coordinates,to=ras, bandwidth=30000)

Error:

UndefVarError: `kerneldensity` not defined

Stacktrace:
 [1] getproperty(x::Module, f::Symbol)
   @ Base ./Base.jl:31
 [2] top-level scope
   @ In[18]:1`

@rafaqz
Copy link
Owner Author

rafaqz commented Dec 24, 2023 •
edited
Loading

Maybe restart your session? If kerneldensity is not defined then you probably don't have this branch active. (if you added that in a session with Rasters already active the new extension code wont be picked up)

Then you also need to do using KernelDensity as wel as using Rasters

The extension part means it only loads when KernelDensity.jl is also loaded.

First maybe try running the code from the tests so we are on the same page:

using Test, Rasters, KernelDensity, Extents

points = collect(zip(rand(100) .* 10, rand(100)))
# Using just res, detecting extent from points
res_rast = kerneldensity(points; res=(0.1, 0.01))
# This doesn't exactly match, someting is wrong with `_extent2dims`
@test_broken map(step, dims(res_rast)) == (0.1, 0.01)
# Using just size, detecting extent from points
size_rast = kerneldensity(points; size=250)
@test size(size_rast) == (250, 250)
# Using a defined extent with res
to = Extents.Extent(X=(0, 10), Y=(0, 5))
ext_res_rast = kerneldensity(points; to, res=0.05)
@test size(ext_rast) == (200, 100)
@test map(step, dims(ext_rast)) == (0.05, 0.05)
# Using a defined extent with size
to = Extents.Extent(X=(0, 10), Y=(0, 5))
ext_res_rast = kerneldensity(points; to, size=(200, 100))
@test size(ext_rast) == (200, 100)
@test map(step, dims(ext_rast)) == (0.05, 0.05)
# Using a defined DimArray
to = rand(X(0:0.01:10), Y(0:0.001:1))
template_rast = kerneldensity(points; to)
@test size(to) == size(template_rast)
@test map(step, dims(template_rast)) == (0.01, 0.001)

# Test plots
using GLMakie
p = Makie.heatmap(res_rast)
Makie.scatter!(points)

p = Makie.heatmap(size_rast)
Makie.scatter!(points)

p = Makie.heatmap(ext_rast)
Makie.scatter!(points)

p = Makie.heatmap(template_rast)
Makie.scatter!(points)

@anjelinejeline
Copy link

Hi @rafaqz
I just tried to run the test but something is missing ( ext_rast pheraphs?)

using Test, Rasters, KernelDensity, Extents

points = collect(zip(rand(100) .* 10, rand(100)))
# Using just res, detecting extent from points
res_rast = kerneldensity(points; res=(0.1, 0.01))
# This doesn't exactly match, someting is wrong with `_extent2dims`
@test_broken map(step, dims(res_rast)) == (0.1, 0.01)
# Using just size, detecting extent from points
size_rast = kerneldensity(points; size=250)
@test size(size_rast) == (250, 250)
# Using a defined extent with res
to = Extents.Extent(X=(0, 10), Y=(0, 5))
ext_res_rast = kerneldensity(points; to, res=0.05)
@test size(ext_rast) == (200, 100)
@test map(step, dims(ext_rast)) == (0.05, 0.05)
# Using a defined extent with size
to = Extents.Extent(X=(0, 10), Y=(0, 5))
ext_res_rast = kerneldensity(points; to, size=(200, 100))
@test size(ext_rast) == (200, 100)
@test map(step, dims(ext_rast)) == (0.05, 0.05)
# Using a defined DimArray
to = rand(X(0:0.01:10), Y(0:0.001:1))
template_rast = kerneldensity(points; to)
@test size(to) == size(template_rast)
@test map(step, dims(template_rast)) == (0.01, 0.001)
Error During Test at In[2]:14
  Test threw exception
  Expression: size(ext_rast) == (200, 100)
  UndefVarError: `ext_rast` not defined
  Stacktrace:
   [1] top-level scope
     @ /software/22.04/julia/julia-1.9.3/share/julia/stdlib/v1.9/Test/src/Test.jl:478
   [2] eval
     @ ./boot.jl:370 [inlined]
   [3] include_string(mapexpr::typeof(REPL.softscope), mod::Module, code::String, filename::String)
     @ Base ./loading.jl:1903
   [4] softscope_include_string(m::Module, code::String, filename::String)
     @ SoftGlobalScope ~/.julia/packages/SoftGlobalScope/u4UzH/src/SoftGlobalScope.jl:65
   [5] execute_request(socket::ZMQ.Socket, msg::IJulia.Msg)
     @ IJulia ~/.julia/packages/IJulia/Vo51o/src/execute_request.jl:67
   [6] #invokelatest#2
     @ ./essentials.jl:819 [inlined]
   [7] invokelatest
     @ ./essentials.jl:816 [inlined]
   [8] eventloop(socket::ZMQ.Socket)
     @ IJulia ~/.julia/packages/IJulia/Vo51o/src/eventloop.jl:8
   [9] (::IJulia.var"#15#18")()
     @ IJulia ./task.jl:514

@rafaqz
Copy link
Owner Author

rafaqz commented Dec 25, 2023

Oh yeah the raster name changed.

@rafaqz
Copy link
Owner Author

rafaqz commented Dec 26, 2023 •
edited
Loading

These are the new fixed tests. You can plot by running specific lines and then the plot commands at the bottom.

using Test, Rasters, KernelDensity, Extents

points = append!(collect(zip(rand(100) .* 10, rand(100))), [(0.0, 0.0), (10.0, 1.0)])

@testset "Using just res, detecting extent from points" begin
    rast = kerneldensity(points; res=(0.1, 0.01))
    @test map(step, dims(rast)) == (0.1, 0.01)
end

@testset "Using just res, detecting extent from points" begin
    rast = kerneldensity(points; size=250)
    @test size(rast) == (250, 250)
end

@testset "Using a defined extent with res" begin
    to = Extents.Extent(X=(0, 10), Y=(0, 5))
    rast = kerneldensity(points; to, res=0.05)
    @test size(rast) == (200, 100)
    @test map(step, dims(rast)) == (0.05, 0.05)
end

@testset "Using a defined extent with size" begin
    to = Extents.Extent(X=(0, 10), Y=(0, 5))
    rast = kerneldensity(points; to, size=(200, 100))
    @test size(rast) == (200, 100)
    @test map(step, dims(rast)) == (0.05, 0.05)
end

@testset "Using a preddefined DimArray" begin
    to = rand(X(0:0.01:10), Y(0:0.001:1))
    rast = kerneldensity(points; to)
    @test size(to) == size(rast)
    @test map(step, dims(rast)) == (0.01, 0.001)
end

# Test plots
# using GLMakie
# p = Makie.heatmap(rast)
# Makie.scatter!(points)

@anjelinejeline
Copy link

Hi @rafaqz I have trouble plotting the test as you suggested

 Info: Precompiling GLMakie [e9467ef8-e4e7-5192-8a1a-b1aee30e663a]
Errors encountered while load FileIO.File{FileIO.DataFormat{:PNG}, String}("/home/panelan/.julia/packages/Makie/Qvk4f/assets/icons/icon-128.png").
All errors:
===========================================
Evaluation into the closed module `ImageIO` breaks incremental compilation because the side effects will not be permanent. This is likely due to some other module mutating `ImageIO` with `eval` during precompilation - don't do this.
===========================================
ArgumentError: Package ImageMagick [6218d12a-5da1-5696-b52f-db25d2ecc6d1] is required but does not seem to be installed:
 - Run `Pkg.instantiate()` to install all recorded dependencies.

===========================================

Fatal error:
ERROR: LoadError: Evaluation into the closed module `ImageIO` breaks incremental compilation because the side effects will not be permanent. This is likely due to some other module mutating `ImageIO` with `eval` during precompilation - don't do this.
Stacktrace:
  [1] eval
    @ ./boot.jl:370 [inlined]
  [2] (::ImageIO.var"#1#2"{Symbol})()
    @ ImageIO ~/.julia/packages/ImageIO/bBdyA/src/ImageIO.jl:19
  [3] lock(f::ImageIO.var"#1#2"{Symbol}, l::ReentrantLock)
    @ Base ./lock.jl:229
  [4] checked_import
    @ ~/.julia/packages/ImageIO/bBdyA/src/ImageIO.jl:8 [inlined]
  [5] load(f::FileIO.File{FileIO.DataFormat{:PNG}, String}; kwargs::Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}})
    @ ImageIO ~/.julia/packages/ImageIO/bBdyA/src/ImageIO.jl:27
  [6] load(f::FileIO.File{FileIO.DataFormat{:PNG}, String})
    @ ImageIO ~/.julia/packages/ImageIO/bBdyA/src/ImageIO.jl:26
  [7] #invokelatest#2
    @ ./essentials.jl:819 [inlined]
  [8] invokelatest
    @ ./essentials.jl:816 [inlined]
  [9] action(::Symbol, ::Vector{Union{Base.PkgId, Module}}, ::FileIO.Formatted; options::Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}})
    @ FileIO ~/.julia/packages/FileIO/BE7iZ/src/loadsave.jl:219
 [10] action
    @ ~/.julia/packages/FileIO/BE7iZ/src/loadsave.jl:196 [inlined]
 [11] action(::Symbol, ::Vector{Union{Base.PkgId, Module}}, ::Symbol, ::String; options::Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}})
    @ FileIO ~/.julia/packages/FileIO/BE7iZ/src/loadsave.jl:185
 [12] action
    @ ~/.julia/packages/FileIO/BE7iZ/src/loadsave.jl:185 [inlined]
 [13] load(::String; options::Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}})
    @ FileIO ~/.julia/packages/FileIO/BE7iZ/src/loadsave.jl:113
 [14] load
    @ ~/.julia/packages/FileIO/BE7iZ/src/loadsave.jl:109 [inlined]
 [15] _broadcast_getindex_evalf
    @ ./broadcast.jl:683 [inlined]
 [16] _broadcast_getindex
    @ ./broadcast.jl:656 [inlined]
 [17] getindex
    @ ./broadcast.jl:610 [inlined]
 [18] copy
    @ ./broadcast.jl:912 [inlined]
 [19] materialize
    @ ./broadcast.jl:873 [inlined]
 [20] icon()
    @ Makie ~/.julia/packages/Makie/Qvk4f/src/Makie.jl:313
 [21] empty_screen(debugging::Bool; reuse::Bool)
    @ GLMakie ~/.julia/packages/GLMakie/WlUOp/src/screen.jl:260
 [22] empty_screen
    @ ~/.julia/packages/GLMakie/WlUOp/src/screen.jl:222 [inlined]
 [23] singleton_screen(debugging::Bool)
    @ GLMakie ~/.julia/packages/GLMakie/WlUOp/src/screen.jl:329
 [24] macro expansion
    @ ~/.julia/packages/GLMakie/WlUOp/src/precompiles.jl:21 [inlined]
 [25] macro expansion
    @ ~/.julia/packages/PrecompileTools/kmH5L/src/workloads.jl:78 [inlined]
 [26] macro expansion
    @ ~/.julia/packages/GLMakie/WlUOp/src/precompiles.jl:18 [inlined]
 [27] macro expansion
    @ ~/.julia/packages/PrecompileTools/kmH5L/src/workloads.jl:140 [inlined]
 [28] top-level scope
    @ ~/.julia/packages/GLMakie/WlUOp/src/precompiles.jl:16
Stacktrace:
  [1] handle_error(e::ErrorException, q::Base.PkgId, bt::Vector{Union{Ptr{Nothing}, Base.InterpreterIP}})
    @ FileIO ~/.julia/packages/FileIO/BE7iZ/src/error_handling.jl:61
  [2] handle_exceptions(exceptions::Vector{Tuple{Any, Union{Base.PkgId, Module}, Vector}}, action::String)
    @ FileIO ~/.julia/packages/FileIO/BE7iZ/src/error_handling.jl:56
  [3] action(::Symbol, ::Vector{Union{Base.PkgId, Module}}, ::FileIO.Formatted; options::Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}})
    @ FileIO ~/.julia/packages/FileIO/BE7iZ/src/loadsave.jl:228
  [4] action
    @ ~/.julia/packages/FileIO/BE7iZ/src/loadsave.jl:196 [inlined]
  [5] action(::Symbol, ::Vector{Union{Base.PkgId, Module}}, ::Symbol, ::String; options::Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}})
    @ FileIO ~/.julia/packages/FileIO/BE7iZ/src/loadsave.jl:185
  [6] action
    @ ~/.julia/packages/FileIO/BE7iZ/src/loadsave.jl:185 [inlined]
  [7] load(::String; options::Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}})
    @ FileIO ~/.julia/packages/FileIO/BE7iZ/src/loadsave.jl:113
  [8] load
    @ ~/.julia/packages/FileIO/BE7iZ/src/loadsave.jl:109 [inlined]
  [9] _broadcast_getindex_evalf
    @ ./broadcast.jl:683 [inlined]
 [10] _broadcast_getindex
    @ ./broadcast.jl:656 [inlined]
 [11] getindex
    @ ./broadcast.jl:610 [inlined]
 [12] copy
    @ ./broadcast.jl:912 [inlined]
 [13] materialize
    @ ./broadcast.jl:873 [inlined]
 [14] icon()
    @ Makie ~/.julia/packages/Makie/Qvk4f/src/Makie.jl:313
 [15] empty_screen(debugging::Bool; reuse::Bool)
    @ GLMakie ~/.julia/packages/GLMakie/WlUOp/src/screen.jl:260
 [16] empty_screen
    @ ~/.julia/packages/GLMakie/WlUOp/src/screen.jl:222 [inlined]
 [17] singleton_screen(debugging::Bool)
    @ GLMakie ~/.julia/packages/GLMakie/WlUOp/src/screen.jl:329
 [18] macro expansion
    @ ~/.julia/packages/GLMakie/WlUOp/src/precompiles.jl:21 [inlined]
 [19] macro expansion
    @ ~/.julia/packages/PrecompileTools/kmH5L/src/workloads.jl:78 [inlined]
 [20] macro expansion
    @ ~/.julia/packages/GLMakie/WlUOp/src/precompiles.jl:18 [inlined]
 [21] macro expansion
    @ ~/.julia/packages/PrecompileTools/kmH5L/src/workloads.jl:140 [inlined]
 [22] top-level scope
    @ ~/.julia/packages/GLMakie/WlUOp/src/precompiles.jl:16
 [23] include(mod::Module, _path::String)
    @ Base ./Base.jl:457
 [24] include(x::String)
    @ GLMakie ~/.julia/packages/GLMakie/WlUOp/src/GLMakie.jl:1
 [25] top-level scope
    @ ~/.julia/packages/GLMakie/WlUOp/src/GLMakie.jl:79
 [26] include
    @ ./Base.jl:457 [inlined]
 [27] include_package_for_output(pkg::Base.PkgId, input::String, depot_path::Vector{String}, dl_load_path::Vector{String}, load_path::Vector{String}, concrete_deps::Vector{Pair{Base.PkgId, UInt128}}, source::String)
    @ Base ./loading.jl:2049
 [28] top-level scope
    @ stdin:3
in expression starting at /home/panelan/.julia/packages/GLMakie/WlUOp/src/precompiles.jl:15
in expression starting at /home/panelan/.julia/packages/GLMakie/WlUOp/src/GLMakie.jl:1
in expression starting at stdin:3
ERROR: LoadError: Failed to precompile GLMakie [e9467ef8-e4e7-5192-8a1a-b1aee30e663a] to "/home/panelan/.julia/compiled/v1.9/GLMakie/jl_6zJ265".
Stacktrace:
  [1] error(s::String)
    @ Base ./error.jl:35
  [2] compilecache(pkg::Base.PkgId, path::String, internal_stderr::IO, internal_stdout::IO, keep_loaded_modules::Bool)
    @ Base ./loading.jl:2300
  [3] compilecache
    @ ./loading.jl:2167 [inlined]
  [4] _require(pkg::Base.PkgId, env::String)
    @ Base ./loading.jl:1805
  [5] _require_prelocked(uuidkey::Base.PkgId, env::String)
    @ Base ./loading.jl:1660
  [6] macro expansion
    @ ./loading.jl:1648 [inlined]
  [7] macro expansion
    @ ./lock.jl:267 [inlined]
  [8] require(into::Module, mod::Symbol)
    @ Base ./loading.jl:1611
  [9] include(fname::String)
    @ Base.MainInclude ./client.jl:478
 [10] top-level scope
    @ REPL[31]:1

@rafaqz
Copy link
Owner Author

rafaqz commented Dec 29, 2023

Did you install ImageMagick like the message says? Do that and try again (there are a few possible image backends so it doesn't install automatically, unfortunately)

@anjelinejeline
Copy link

Yes, thanks this issue can be fixed if ImageMagick is installed .. but I am now getting an error related to rast which is not defined apparently

@rafaqz
Copy link
Owner Author

rafaqz commented Dec 31, 2023 •
edited
Loading

You need to run the specific line individually not in the test block just copy and paste it out

The run the plot code for whatever youre interested in

@anjelinejeline
Copy link

Hi rafa .. it works, but I cannot plot it .. this code returns a Scatter{Tuple{Vector{Point{2, Float32}}}}

using GLMakie, Makie
p = Makie.heatmap(rast)
Makie.scatter!(points)

@rafaqz
Copy link
Owner Author

rafaqz commented Jan 2, 2024 •
edited
Loading

Not sure how you are running these things? It should call display automatically if you do that in the REPL.

But you can call display manually on the output if thats what you want. display(p)

(Are you always running things as whole scripts? the best way to work with Julia is in a live REPL session)

@anjelinejeline
Copy link

anjelinejeline commented Jan 3, 2024 •
edited
Loading

Hi rafa .. it works!! thanks .. the output of the function also seems OKAY

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Reviewers
No reviews
Assignees
No one assigned
Labels
None yet
Projects
None yet
Milestone
No milestone
Development

Successfully merging this pull request may close these issues.
2 participants
@rafaqz @anjelinejeline