FCUP

• FCUP
• Define the TMC
• Define the seeds
• Compute the streamlines
• plot the streamlines
• Compute the connections

This is a more advanced tutorial because we want to show how to apply a mask.

Define the TMC

import Tractography as TG

model = TG.TMC(Δt = 0.125f0,
            foddata = TG.FODData((@__DIR__) * "/../../examples/fod-FC.nii.gz"),
            cone = TG.Cone(45),
            proba_min = 0.015f0,
            )

TMC with elype Float32
 ├─ Δt = 0.125
 ├─ minimal probability = 0.015
 ├─ cone                = Cone{Int64}(45)
 ├─ mollifier           = max_mollifier
 ├─ evaluation of SH    = PreComputeAllFOD()
 └─ data : (lmax = 8)

Just for fun, we plot the FODs of the model.

using CairoMakie

f, sc = TG.plot_fod(model; n_sphere = 1500, radius = 0.3, st = 2);
cam3d = Makie.cameracontrols(sc)
cam3d.eyeposition[] = Vec3f(85, 95, -28)
cam3d.lookat[] = Vec3f(84, 95, 59)
rotate_cam!(sc.scene, 0, 0, -pi/2)
f

Define the seeds

We next apply a mask on the boundary of which the streamlines stop.

using NIfTI
mask = NIfTI.niread((@__DIR__) * "/../../examples/wm-FC.nii.gz");
TG._apply_mask!(model, mask);

┌ Warning: #= line 0 =#:
│ `LoopVectorization.check_args` on your inputs failed; running fallback `@inbounds @fastmath` loop instead.
│ Use `warn_check_args=false`, e.g. `@turbo warn_check_args=false ...`, to disable this warning.
└ @ LoopVectorization ~/.julia/packages/LoopVectorization/GKxH5/src/condense_loopset.jl:1166

We compute Nmc streamlines, hence we need Nmc seeds

Nmc = 100_000
seeds = TG.from_fod(model, Nmc; maxfod_start = true)

6×100000 Matrix{Float32}:
 36.3744    84.2544    35.9498     81.4211    …  105.494       120.408
 69.3238    63.0655    68.8369     33.1784        62.8183      117.472
  6.48182    3.12279    8.74504     8.63106        6.01086       2.53544
 -0.959548   0.7217    -0.997282    0.629208      -0.661177      0.73366
  0.165176  -0.684349  -0.0653648   0.776403      -0.750227     -0.113328
 -0.228      0.104      0.034      -0.036     …   -0.00200002    0.67

Compute the streamlines

streamlines, tract_length = TG.sample(model, TG.Deterministic(), seeds; nt = 1000);
println("Dimension of computed streamlines = ", size(streamlines))

kernel : 6.608052 seconds (984.76 k allocations: 48.313 MiB, 7.38% compilation time)
Dimension of computed streamlines = (5, 1000, 100000)

plot the streamlines

f, scene = @time TG.plot_fod(model; n_sphere = 500, radius = 0.3, st = 1);
ind_st = findall(tract_length .> 60)
TG.plot_streamlines!(scene, streamlines[:, :, ind_st[1:10:end]])
f

We can also add the seeds

scatter!(scene, seeds[1:3, ind_st[1:10:end]], color = :white)
f

Compute the connections

When computing structural connectivity, we don't need to record the entire streamline but only its extremities.

streamlines, tract_length = TG.sample(model, TG.Connectivity(TG.Deterministic()), seeds; nt = 1000);
println("Dimension of computed streamlines = ", size(streamlines))

kernel : 6.064950 seconds (405.28 k allocations: 19.390 MiB, 6.40% compilation time)
Dimension of computed streamlines = (5, 2, 100000)