Questions about CD_ELECTRODE judgement

[zql2021]

Hi Felix,
I’m encountering an issue where the end position of charge are incorrectly flagged as CD_ELECTRODE during simulating a gap-event in a double-sided strip HPGe detector.

I added println in function _check_and_update_position!to debug

The output shows：

But the charge actually doesn‘t’ reach the electrode (x-direction)：

the simulation code:

T=Float64
position = (0.0, 1.5, 2.5)  # mm
pos = CartesianPoint{T}(position[1] / 1000.0, position[2] / 1000.0, position[3] / 1000.0)# convert mm to m
Edep = T(100) * u"keV"
evt=SolidStateDetectors.Event(pos,Edep)
simulate!(evt, sim, max_nsteps = 10000, Δt =0.1u"ns", diffusion=true,self_repulsion=true,verbose = true)
plot(
    plot( evt.waveforms[5] ,lw = 1.5, xlims = (0, 200), xlabel = "Time", ylabel = "Charge", 
unitformat = :slash, legend = false, tickfontsize = 12, guidefontsize = 14,title="ac_leftimage"),
    plot( evt.waveforms[6] ,lw = 1.5, xlims = (0, 200), xlabel = "Time", ylabel = "Charge",
unitformat = :slash, legend = false, tickfontsize = 12, guidefontsize = 14,title="ac_collect"),
    plot( evt.waveforms[7] ,lw = 1.5, xlims = (0, 200), xlabel = "Time", ylabel = "Charge",
unitformat = :slash, legend = false, tickfontsize = 12, guidefontsize = 14,title="ac_rightimage"),
    plot( evt.waveforms[14] ,lw = 1.5, xlims = (0, 200), xlabel = "Time", ylabel = "Charge",
unitformat = :slash, legend = false, tickfontsize = 12, guidefontsize = 14,title="dc_leftimage"),
    plot( evt.waveforms[15] ,lw = 1.5, xlims = (0, 200), xlabel = "Time", ylabel = "Charge",
unitformat = :slash, legend = false, tickfontsize = 12, guidefontsize = 14,title="dc_collect"),
    plot( evt.waveforms[16] ,lw = 1.5, xlims = (0, 200), xlabel = "Time", ylabel = "Charge",
unitformat = :slash, legend = false, tickfontsize = 12, guidefontsize = 14,title="dc_rightimage"),
size = (1400, 800), layout = (2, 3)
)

The waveform shows incomplete charge collection on the ac_collect electrode (because weighting potential of the end position is less than 1), while neighboring electrodes collect partial ：

Further check of function get_crossing_pos：

for contact in det.contacts  
    for surf in ConstructiveSolidGeometry.surfaces(contact.geometry)  
        for pt in ConstructiveSolidGeometry.intersection(surf, line)  
            if pt in contact &&  
                0 ≤ (pt - pt_in) ⋅ direction ≤ 1 &&  # pt between pt_in and pt_out  
                norm(pt - pt_in) < norm(crossing_pos[1] - pt_in)  
                crossing_pos = (pt, CD_ELECTRODE, normalize(ConstructiveSolidGeometry.normal(surf, pt)))  
                println("crossing_pos: ", crossing_pos)  
            end  
        end  
    end  
end  

I found pt_in and pt_out, their x-coordinates don't reach the electrode edge:

For the issue described above, could you give me some suggestions to fix the termination checks and let the charge end at the electrode edge?
Below is the detector config file (detector.yaml) for reference (if possble, could you check the contact geometry? I didn't find any issues when I check).
detector.zip
Looking forward to your reply, thank you very much！

[fhagemann]

How did you calculate the electric potential, field and weighting potential?
Using the standard simulate!? What refinement did you use? I will try to reproduce this and to look into this, once I find the time.

[zql2021]

I calculate the electric field ,potential by the code below

 T = Float64
    sim = Simulation{T}(detconfigPath)
    simulate!(sim, min_tick_distance=0.001u"mm", verbose=true)
    calculate_electric_potential!(sim, 
        refinement_limits=[0.2, 0.1, 0.05, 0.02], 
    )
    calculate_electric_field!(sim)

[fhagemann]

This issue could result from the tolerance that we apply to intersections, see L. 346 here:

SolidStateDetectors.jl/src/ChargeDrift/ChargeDrift.jl

Lines 319 to 368 in 5f3a5cb

	function get_crossing_pos( det::SolidStateDetector{T}, point_types::PointTypes{T, 3, S}, pt_in::CartesianPoint{T}, pt_out::CartesianPoint{T};
	max_n_iter::Int = 500)::Tuple{CartesianPoint{T}, UInt8, CartesianVector{T}} where {T <: SSDFloat, S}
	# check if the points are already in contacts
	if pt_in in det.contacts return (pt_in, CD_ELECTRODE, CartesianVector{T}(0,0,0)) end
	direction::CartesianVector{T} = normalize(pt_out - pt_in)
	crossing_pos::Tuple{CartesianPoint{T}, UInt8, CartesianVector{T}} = (pt_out, CD_OUTSIDE, CartesianVector{T}(0,0,0)) # need undef version for this
	# define a Line between pt_in and pt_out
	line::ConstructiveSolidGeometry.Line{T} = ConstructiveSolidGeometry.Line{T}(pt_in, direction)
	# check if the Line intersects with a surface of a Contact
	for contact in det.contacts
	for surf in ConstructiveSolidGeometry.surfaces(contact.geometry)
	for pt in ConstructiveSolidGeometry.intersection(surf, line)
	if pt in contact &&
	0 ≤ (pt - pt_in) ⋅ direction ≤ 1 && # pt within pt_in and pt_out
	norm(pt - pt_in) < norm(crossing_pos[1] - pt_in) # pt closer to pt_in that previous crossing_pos
	crossing_pos = (pt, CD_ELECTRODE, normalize(ConstructiveSolidGeometry.normal(surf, pt)))
	end
	end
	end
	end
	# if the Line does not intersect with a surface of a Contact, check if it does intersect with the surface of the Semiconductor
	if crossing_pos[2] & CD_OUTSIDE > 0
	tol::T = 5000 * ConstructiveSolidGeometry.csg_default_tol(T)
	# check if the Line intersects with a surface of the Semiconductor
	for surf in ConstructiveSolidGeometry.surfaces(det.semiconductor.geometry)
	for pt in ConstructiveSolidGeometry.intersection(surf, line)
	normal = normalize(ConstructiveSolidGeometry.normal(surf, pt))
	if pt + tol * normal in det.semiconductor && # point "before" crossing_pos should be in
	!(pt - tol * normal in det.semiconductor) && # point "after" crossing_pos should be out
	0 ≤ (pt - pt_in) ⋅ direction ≤ 1 && # pt within pt_in and pt_out
	norm(pt - pt_in) < norm(crossing_pos[1] - pt_in) # pt closer to pt_in that previous crossing_pos
	CD_POINTTYPE::UInt8 = point_types[pt] & update_bit == 0 ? CD_ELECTRODE : CD_FLOATING_BOUNDARY
	crossing_pos = (pt, CD_POINTTYPE, normal)
	end
	end
	end
	end
	# if there is no intersection, check if the next point is in (within the tolerance)
	if (crossing_pos[2] & CD_OUTSIDE > 0) && pt_out in det.contacts
	return (pt_out, CD_ELECTRODE, CartesianVector{T}(0,0,0))
	end
	crossing_pos
	end

Here, we chose a multiplication factor of 5000 with the default tolerance to account for rounding issues.
You could decrease this factor or remove it completely to see whether that fixes your problem.
However, this might cause other (unexpected) problems, so I would be careful removing it from SSD in general.

Also, it thinks it is in the electrode because we check for the point_types here.
To avoid costly geometry checks, we first check if a point has a fixed point type, meaning that it does not change its value in the potential calculation because it has a fixed potential applied to it.
The point types might over-extend the contact into the bulk.

In your case with your simulation settings:
You can check the grid points that your x axis has (in the region 0 mm < x < 5 mm):

filter(t -> 0 < t < 0.5, (sim.point_types.grid.x.ticks)*1000)

4-element Vector{Float64}:
 0.0625000000000007
 0.1250000000000014
 0.18750000000000072
 0.25

The center between grid points can be calculated with StatsBase.midpoints:

using StatsBase
filter(t -> 0 < t < 0.5, StatsBase.midpoints((sim.point_types.grid.x.ticks)*1000))

4-element Vector{Float64}:
 0.03125000000000035
 0.09375000000000104
 0.15625000000000105
 0.21875000000000036

which in your case means that everything with x > 0.21875mm will be flagged as "ELECTRODE" in our charge drift code.

In your case, this is not optimal.. So if you come up with an elegant solution that does not break any of the current tests, we'd be happy to incorporate that in the software :)

[zql2021]

Hi felix ,Thank you so much for your detailed response during the holiday, I will try the approaches you suggested.