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}:
 116.694      74.7159    120.024     …  45.0422     36.0795     141.105
  26.7132    117.01      119.678        71.6846     53.9805     119.771
   8.87195     6.17824     9.15545       6.19173     8.8967       2.93544
  -0.689921   -0.405147   -0.977391     -0.994637    0.997395    -0.982866
   0.710609    0.906075   -0.197308      0.0171193   0.0240615   -0.11506
   0.138       0.122       0.076     …   0.102       0.068        0.144

Compute the streamlines

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

kernel : 7.229410 seconds (961.05 k allocations: 47.428 MiB, 6.58% compilation time)
Dimension of computed streamlines = (3, 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.138498 seconds (406.88 k allocations: 19.530 MiB, 6.16% compilation time)
Dimension of computed streamlines = (3, 2, 100000)