diff --git a/Project.toml b/Project.toml
index 7e3b2f3..8e15430 100644
--- a/Project.toml
+++ b/Project.toml
@@ -10,11 +10,15 @@ Gtk4 = "9db2cae5-386f-4011-9d63-a5602296539b"
 GtkObservables = "8710efd8-4ad6-11eb-33ea-2d5ceb25a41c"
 ImageCore = "a09fc81d-aa75-5fe9-8630-4744c3626534"
 ImageIO = "82e4d734-157c-48bb-816b-45c225c6df19"
+ImageMagick = "6218d12a-5da1-5696-b52f-db25d2ecc6d1"
 ImageMorphology = "787d08f9-d448-5407-9aad-5290dd7ab264"
 ImageSegmentation = "80713f31-8817-5129-9cf8-209ff8fb23e1"
 ImageView = "86fae568-95e7-573e-a6b2-d8a6b900c9ef"
 JLD2 = "033835bb-8acc-5ee8-8aae-3f567f8a3819"
+LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
+OpenCV = "f878e3a2-a245-4720-8660-60795d644f2a"
 Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c"
+Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2"
 XLSX = "fdbf4ff8-1666-58a4-91e7-1b58723a45e0"
 
 [compat]
@@ -24,10 +28,12 @@ Gtk4 = "0.7"
 GtkObservables = "2"
 ImageCore = "0.10"
 ImageIO = "0.6"
+ImageMagick = "1.4.2"
 ImageMorphology = "0.4"
 ImageSegmentation = "1.9"
 ImageView = "0.12"
 JLD2 = "0.5"
+OpenCV = "4.6.1"
 XLSX = "0.10"
 julia = "1.10"
 
diff --git a/docs/src/assets/Picture.png b/docs/src/assets/Picture.png
index bf9f5c5..a0e56dd 100644
Binary files a/docs/src/assets/Picture.png and b/docs/src/assets/Picture.png differ
diff --git a/docs/src/assets/blurred_calibration.png b/docs/src/assets/blurred_calibration.png
new file mode 100644
index 0000000..f0760d9
Binary files /dev/null and b/docs/src/assets/blurred_calibration.png differ
diff --git a/docs/src/assets/gui.png b/docs/src/assets/gui.png
index 3256fd2..58e2658 100644
Binary files a/docs/src/assets/gui.png and b/docs/src/assets/gui.png differ
diff --git a/docs/src/index.md b/docs/src/index.md
index 54ff3fc..c156c29 100644
--- a/docs/src/index.md
+++ b/docs/src/index.md
@@ -13,7 +13,7 @@ The yellow spot corresponds to a stimulus provided by the experimenter, and the
 
 Tips on image quality:
 
-- Put the stimulus near one of the four corners
+- Put the stimulus near one of the four corners, keeping its location as consistent as possible between images
 - Ensure lighting is fairly uniform
 - Make sure that any extraneous marks (e.g., the black writing in the image above) are of a very different color from scent marks.
 - Ensure that all your images are of the same size (i.e., same number of pixels horizontally and vertically), even if there are some extra pixels on the edges of the image
@@ -99,12 +99,37 @@ Alternatively, you can supply a list of files:
 julia> gui("results_file_name", ["PictureA.png", "mouse7.png"])
 ```
 
+Additionally, you can supply a calibration image to improve segmentation by correcting for uneven illumination.
+
+![An image to be used for calibration](assets/blurred_calibration.png)
+
+```
+julia> gui("results_file_name", ["PictureA.png", "mouse7.png"]; background_path="calibration_image.png")
+```
+
+There are a few more keyword arguments that might be useful:
+  - `crop_top`, `crop_bottom`, `crop_left`, and `crop_right` specify a number of pixels to be cropped from images along each side
+  - `expectedloc` specifies the expected location of the stimulus spot in pixels (after cropping)
+
+```
+julia> gui("results_file_name", glob"Picture*.png"; 
+            background_path="calibration_image.png,
+            crop_top=93,
+            crop_bottom=107,
+            crop_left=55,
+            crop_right=45,    
+            expectedloc=[1600,3333]
+        )
+```
+
 However you launch it, you should see something like this:
 
 ![GUI](assets/gui.png)
 
 On the top is the raw image. On the bottom is the segmented image; you should visually compare the two to check whether you're pleased with the quality of the segmentation.
-(If not, click "Skip" to omit that file from analysis.)
+
+If the default segmentation doesn't look quite right, try adjusting the Color Similarity Threshold value using the buttons at the bottom of the GUI.
+(If you can't obtain a segmentation that you're happy with, click "Skip" to omit that file from analysis.)
 
 If you like the segmentation, your tasks are:
 - click on all the checkboxes with colors that correspond to urine spots. You'll notice that the stimulus spot is pre-clicked (you can correct its choice if it didn't pick correctly). Most of the time there will be only one you need to check, but you can click more than one.
@@ -124,6 +149,8 @@ julia> using CounterMarking, ImageView   # load packages (if this is a fresh ses
 
 julia> count = density_map("results_file_name.jld2");
 
+julia> dmap = count ./ maximum(count);
+
 julia> dct = imshow(dmap);
 ```
 
diff --git a/src/CounterMarking.jl b/src/CounterMarking.jl
index edd9eff..93d7b7e 100644
--- a/src/CounterMarking.jl
+++ b/src/CounterMarking.jl
@@ -11,6 +11,8 @@ using Gtk4
 using GtkObservables
 using ImageView
 using Random
+using LinearAlgebra
+using Statistics
 
 export segment_image, stimulus_index, spots, Spot, upperleft
 export writexlsx, process_images, density_map
diff --git a/src/gui.jl b/src/gui.jl
index 32d76a7..5700f08 100644
--- a/src/gui.jl
+++ b/src/gui.jl
@@ -20,15 +20,31 @@ The JLD2 file can be used by [`density_map`](@ref).
 """
 function gui(
         outbase::AbstractString, files;
-        colors=distinguishable_colors(15, [RGB(1, 1, 1)]; dropseed=true),
+        colors=distinguishable_colors(15, [RGB(1, 1, 1)]; dropseed=false),
         btnclick = Condition(),         # used for testing
         whichbutton = Ref{Symbol}(),    # used for testing
         preclick::Union{Int,Nothing} = nothing,  # used for testing
+        background_path = nothing,  # used to correct for non-uniform illumination (see joinpath(pkgdir(@__MODULE__),"docs","src","assets","blurred_calibration.png"))
+        crop_top::Int = 0,          # crop this many pixels off of each side
+        crop_bottom::Int = 0,
+        crop_left::Int = 0,
+        crop_right::Int = 0,
+        colorproj = RGB{Float32}(1, 1, -2), # used for identifying the stimulus
+        expectedloc = nothing,              # ""
     )
     channelpct(x) = string(round(Int, x * 100)) * '%'
 
     outbase, _ = splitext(outbase)
 
+    # Initialize segmented image and color similarity threshold
+    img = nothing
+    rescaledimg = nothing
+    bkgimg = isnothing(background_path) ? nothing : Float32.(Gray.(load(background_path)[crop_top+1:end-crop_bottom, crop_left+1:end-crop_right]))
+    bkgmean = isnothing(bkgimg) ? nothing : Float32(mean(bkgimg))
+    seg = nothing
+    threshold = 0.15
+    labels2idx = Dict{Int,Int}()
+    idx2labels = Dict{Int,Int}()
     # Set up basic properties of the window
     winsize = round.(Int, 0.8 .* screen_size())
     win = GtkWindow("CounterMarking", winsize...)
@@ -86,6 +102,64 @@ function gui(
         seggrid[col, row] = cb.widget
         push!(cbs, cb)
     end
+
+    thrshbx = GtkBox(:h)
+    thrshlbl = GtkLabel("Color Similarity Threshold: $(threshold)")
+    push!(thrshbx, thrshlbl)
+    thrshbtnbx = GtkBox(:v)
+    push!(thrshbx, thrshbtnbx)
+    lgincbtn = GtkButton("+0.05")
+    smincbtn = GtkButton("+0.01")
+    smdecbtn = GtkButton("-0.01")
+    lgdecbtn = GtkButton("-0.05")
+    push!(thrshbtnbx, lgincbtn)
+    push!(thrshbtnbx, smincbtn)
+    push!(thrshbtnbx, smdecbtn)
+    push!(thrshbtnbx, lgdecbtn)
+    push!(guibx, thrshbx)
+
+    update_threshold = change -> begin
+        newvalue = round(threshold + change; digits=2)
+        try
+            seg = segment_image(rescaledimg; threshold = newvalue)
+            nsegs = length(segment_labels(seg))
+            nsegs > length(colors) && @warn "More than $(length(colors)) segments ($(nsegs)). Excluded ones will be displayed in white and will not be selectable"
+            empty!(labels2idx)
+            empty!(idx2labels)
+            for (i,l) in enumerate(sort(seg.segment_labels; rev=true, by=(l -> segment_pixel_count(seg,l))))
+                idx = i == 1 ? 2 : i == 2 ? 1 : i
+                idx = i > 15 ? 1 : idx  # show remaining segments in white
+                push!(labels2idx, l=>idx)
+                push!(idx2labels, idx=>l)
+            end
+            centroidsacc, _ = get_centroidsacc(seg.image_indexmap)
+            istim = labels2idx[stimulus_index(seg, centroidsacc; colorproj = colorproj, expectedloc = expectedloc)]            
+            for (j, cb) in enumerate(cbs)
+                # set_gtk_property!(cb, "active", j <= nsegs)
+                cb[] = (j == istim || j == preclick)
+            end
+            imshow(canvases[1, 1], img)
+            imshow(canvases[2, 1], map(i->colors[labels2idx[i]], labels_map(seg)))
+            threshold = newvalue
+            thrshlbl.label = "Color Similarity Threshold: $threshold"
+        catch e
+            @show e
+            # display?
+        end
+    end
+    signal_connect(lgincbtn, "clicked") do w, others...
+        update_threshold(0.05)
+    end
+    signal_connect(smincbtn, "clicked") do w, others...
+        update_threshold(0.01)
+    end
+    signal_connect(smdecbtn, "clicked") do w, others...
+        update_threshold(-0.01)
+    end
+    signal_connect(lgdecbtn, "clicked") do w, others...
+        update_threshold(-0.05)
+    end
+
     # Add "Done & Next" and "Skip" buttons
     donebtn = button("Done & Next")
     skipbtn = button("Skip")
@@ -102,17 +176,30 @@ function gui(
 
     results = []
     for (i, file) in enumerate(files)
-        img = color.(load(file))
-        seg = segment_image(img)
+        threshold = 0.15
+        thrshlbl.label = "Color Similarity Threshold: $threshold"
+        img = color.(load(file))[crop_top+1:end-crop_bottom, crop_left+1:end-crop_right]
+        rescaledimg = isnothing(bkgimg) ? img : (img ./ bkgimg .* bkgmean)
+
+        seg = segment_image(rescaledimg; threshold = threshold)
         nsegs = length(segment_labels(seg))
-        @assert nsegs < length(colors) "Too many segments for colors"
-        istim = stimulus_index(seg)
+        nsegs > length(colors) && @warn "More than $(length(colors)) segments ($(nsegs)). Excluded ones will be displayed in white and will not be selectable"
+        empty!(labels2idx)
+        empty!(idx2labels)
+        for (i,l) in enumerate(sort(seg.segment_labels; rev=true, by=(l -> segment_pixel_count(seg,l))))
+            idx = i == 1 ? 2 : i == 2 ? 1 : i
+            idx = i > 15 ? 1 : idx  # show remaining segments in white
+            push!(labels2idx, l=>idx)
+            push!(idx2labels, idx=>l)
+        end
+        centroidsacc, _ = get_centroidsacc(seg.image_indexmap)
+        istim = labels2idx[stimulus_index(seg, centroidsacc; colorproj = colorproj, expectedloc = expectedloc)]            
         for (j, cb) in enumerate(cbs)
             # set_gtk_property!(cb, "active", j <= nsegs)
             cb[] = (j == istim || j == preclick)
-        end
+        end    
         imshow(canvases[1, 1], img)
-        imshow(canvases[2, 1], map(i->colors[i], labels_map(seg)))
+        imshow(canvases[2, 1], map(l->colors[labels2idx[l]], labels_map(seg)))
 
         wait(btnclick)
         whichbutton[] == :skip && continue
@@ -120,13 +207,13 @@ function gui(
         keep = Int[]
         for (j, cb) in enumerate(cbs)
             if cb[]
-                push!(keep, j)
+                push!(keep, idx2labels[j])
             end
         end
         pixelskeep = map(i -> i ∈ keep, labels_map(seg))
         L = label_components(pixelskeep)
         newseg = SegmentedImage(img, L)
-        spotdict, stimulus = spots(newseg)
+        spotdict, stimulus = spots(newseg; colorproj = colorproj, expectedloc = expectedloc)
         push!(results, (file, spotdict, stimulus, newseg))
     end
 
diff --git a/src/segment.jl b/src/segment.jl
index ece05d6..fbdff04 100644
--- a/src/segment.jl
+++ b/src/segment.jl
@@ -8,9 +8,9 @@ Larger `threshold` values will result in fewer segments.
 """
 function segment_image(
         img::AbstractMatrix{<:Color};
-        threshold::Real = 0.2,      # threshold for color similarity in region growing
+        threshold::Real = 0.15,     # threshold for color similarity in region growing
         prune::Bool = true,         # prune small segments
-        min_size::Int = 50,         # minimum size of segments to keep
+        min_size::Int = 500,        # minimum size of segments to keep
     )
     seg = unseeded_region_growing(img, threshold)
     if prune
@@ -22,16 +22,30 @@ end
 segment_image(img::AbstractMatrix{<:Colorant}; kwargs...) = segment_image(color.(img); kwargs...)
 
 """
-    idx = stimulus_index(seg::SegmentedImage, colorproj = RGB(1, 1, -2))
+    idx = stimulus_index(seg::SegmentedImage, centroidsacc; colorproj = RGB(1, 1, -2), expectedloc = nothing)
 
 Given a segmented image `seg`, return the index of the segment that scores
 highest on the product of (1) projection (dot product) with `colorproj` and (2)
 number of pixels.
+
+Optionally, if images were taken with a fixed location for the stimulus, a segment's score
+is divided by the squared distance of its centroid (via `centroidsacc`) from the position given by `expectedloc`.
 """
-function stimulus_index(seg::SegmentedImage, colorproj = RGB{Float32}(1, 1, -2))
-    proj = [l => (colorproj ⋅ segment_mean(seg, l)) * segment_pixel_count(seg, l) for l in segment_labels(seg)]
-    (i, _) = argmax(last, proj)
-    return i
+function stimulus_index(seg::SegmentedImage, centroidsacc; colorproj = RGB{Float32}(1, 1, -2), expectedloc = nothing)
+    if !isnothing(expectedloc)
+        proj = map(segment_labels(seg)) do l
+            l == 0 && return 0
+            val = centroidsacc[l]
+            centroid = [round(Int, val[1] / val[3]), round(Int, val[2] / val[3])]
+            return l => (colorproj ⋅ segment_mean(seg, l) * segment_pixel_count(seg, l) / max(1, sum(abs2, centroid .- expectedloc)))
+        end
+        (i, _) = argmax(last, proj)
+        return i
+    else 
+        proj = [l => (colorproj ⋅ segment_mean(seg, l)) * segment_pixel_count(seg, l) for l in segment_labels(seg)]
+        (i, _) = argmax(last, proj)
+        return i
+    end
 end
 
 # function contiguous(seg::SegmentedImage, img::AbstractMatrix{<:Color}; min_size::Int = 50)
@@ -47,6 +61,35 @@ end
 # contiguous(seg::SegmentedImage, img::AbstractMatrix{<:Colorant}; kwargs...) =
 #     contiguous(seg, color.(img); kwargs...)
 
+"""
+    centroidsacc, nadj = get_centroidsacc(seg::SegmentedImage)
+
+Given a the index map `indexmap` of a segmented image, return an accumulator for each segment's centroid
+as well as the number of times two segments are adjacent.
+"""
+function get_centroidsacc(indexmap::Matrix{Int64})
+    keypair(i, j) = i < j ? (i, j) : (j, i)
+    R = CartesianIndices(indexmap)
+    Ibegin, Iend = extrema(R)
+    I1 = oneunit(Ibegin)
+    centroidsacc = Dict{Int, Tuple{Int, Int, Int}}()   # accumulator for centroids
+    nadj = Dict{Tuple{Int, Int}, Int}()             # number of times two segments are adjacent
+    for idx in R
+        l = indexmap[idx]
+        l == 0 && continue
+        acc = get(centroidsacc, l, (0, 0, 0))
+        centroidsacc[l] = (acc[1] + idx[1], acc[2] + idx[2], acc[3] + 1)
+        for j in max(Ibegin, idx - I1):min(Iend, idx + I1)
+            lj = indexmap[j]
+            if lj != l && lj != 0
+                k = keypair(l, lj)
+                nadj[k] = get(nadj, k, 0) + 1
+            end
+        end
+    end
+    return centroidsacc, nadj
+end
+
 struct Spot
     npixels::Int
     centroid::Tuple{Int, Int}
@@ -67,40 +110,21 @@ Spots larger than `max_size_frac * npixels` (default: 10% of the image) are igno
 function spots(
         seg::SegmentedImage;
         max_size_frac=0.1,            # no spot is bigger than max_size_frac * npixels
+        kwargs...
     )
-    keypair(i, j) = i < j ? (i, j) : (j, i)
-
-    istim = stimulus_index(seg)
+    centroidsacc, nadj = get_centroidsacc(seg.image_indexmap)
+    istim = stimulus_index(seg, centroidsacc; kwargs...)
 
-    label = seg.image_indexmap
-    R = CartesianIndices(label)
-    Ibegin, Iend = extrema(R)
-    I1 = oneunit(Ibegin)
-    centroidsacc = Dict{Int, Tuple{Int, Int, Int}}()   # accumulator for centroids
-    nadj = Dict{Tuple{Int, Int}, Int}()             # number of times two segments are adjacent
-    for idx in R
-        l = label[idx]
-        l == 0 && continue
-        acc = get(centroidsacc, l, (0, 0, 0))
-        centroidsacc[l] = (acc[1] + idx[1], acc[2] + idx[2], acc[3] + 1)
-        for j in max(Ibegin, idx - I1):min(Iend, idx + I1)
-            lj = label[j]
-            if lj != l && lj != 0
-                k = keypair(l, lj)
-                nadj[k] = get(nadj, k, 0) + 1
-            end
-        end
-    end
     stimulus = Ref{Pair{Int,Spot}}()
     filter!(centroidsacc) do (key, val)
         if key == istim
             stimulus[] = key => Spot(val[3], (round(Int, val[1] / val[3]), round(Int, val[2] / val[3])))
             return false
         end
-        val[3] <= max_size_frac * length(label) || return false
+        val[3] <= max_size_frac * length(seg.image_indexmap) || return false
         # # is the centroid within the segment?
         # x, y = round(Int, val[1] / val[3]), round(Int, val[2] / val[3])
-        # l = label[x, y]
+        # l = seg.image_indexmap[x, y]
         # @show l
         # l == key || return false
         # is the segment lighter than most of its neighbors?
diff --git a/test/runtests.jl b/test/runtests.jl
index a8ae579..e5e4044 100644
--- a/test/runtests.jl
+++ b/test/runtests.jl
@@ -1,13 +1,22 @@
 using CounterMarking
 using FileIO
 using XLSX
+using ImageCore
 using ImageView
 using Glob
 using Test
+using Statistics
 
 @testset "CounterMarking.jl" begin
     testdir = joinpath(pkgdir(CounterMarking), "docs", "src", "assets")
-    img = load(joinpath(testdir, "Picture.png"))
+    oimg = load(joinpath(testdir, "Picture.png"))
+
+    # settings particular to the test image:
+    crop_top, crop_bottom, crop_left, crop_right = 93, 107, 55, 45
+    bkgimg = Float32.(Gray.(load(joinpath(testdir,"blurred_calibration.png"))[crop_top+1:end-crop_bottom, crop_left+1:end-crop_right]))
+    expectedloc = [1600,3333]
+
+    img = (oimg[crop_top+1:end-crop_bottom, crop_left+1:end-crop_right] ./ bkgimg .* Float32(mean(bkgimg)))
     seg = segment_image(img)
     dct = meanshow(seg)
     @test haskey(dct, "gui")
@@ -19,11 +28,11 @@ using Test
     @test haskey(dct, "window")
     ImageView.closeall()
 
-    spotdict, stimulus = spots(seg)
+    spotdict, stimulus = spots(seg; expectedloc=expectedloc)
     _, stimspot = stimulus
     @test stimspot.npixels > 1000
-    @test stimspot.centroid[1] < size(img, 1) ÷ 2
-    @test stimspot.centroid[2] > size(img, 2) ÷ 2
+    @test stimspot.centroid[1] - expectedloc[1] < 50
+    @test stimspot.centroid[2] - expectedloc[2] < 50
 
     stdspotdict, stdstimulus = upperleft(spotdict, stimulus, size(img))
     _, stimspot = stdstimulus