Donders-Institute / PRESTUS

PREprocessing & Simulations for Transcranial Ultrasound Stimulation package
GNU General Public License v3.0
15 stars 10 forks source link

Acoustic simulations with NaNs #48

Open jkosciessa opened 1 month ago

jkosciessa commented 1 month ago

With the updated tissue attenuation, I run into sporadic acoustic simulations that are filled with NaNs. Those will eventually crash and are characterized by max ISPPA positions of 1 1 1:

Final transducer, expected focus, and max ISPPA positions

ans =

   175   175     1
   189   189     1
   363   194     1

The crash will happen at the attempt to plot, and will be of the format:

{Error using matlab.graphics.axis.Axes/set
Error setting property 'CLim' of class 'Axes':
Value must be a 1x2 vector of numeric type in which the second element is
greater than the first element or is Inf.

Error in clim (line 142)
        set(axHandle, 'CLim',  climits);

Error in plot_isppa_over_image (line 201)
        clim(options.isppa_color_range);

Error in single_subject_pipeline (line 379)
        [~,~,~,~,~,~,~,h]=plot_isppa_over_image(...

Error in tpb5b9b4be_c0ba_4372_a00e_bf852ef7b5b1 (line 1)
load /project/2424103.01/thalstim_simulations/thalstim_sim/data/tussim/CTX500-026-010_79.6mm_pCT_60W_debug/sub-007/batch_job_logs/tpc8e4bab0_f2a1_40a7_9c57_54ef6096e020.mat; cd /project/2424103.01/thalstim_simulations/thalstim_sim/tools/PRESTUS; single_subject_pipeline(subject_id, parameters); delete /project/2424103.01/thalstim_simulations/thalstim_sim/data/tussim/CTX500-026-010_79.6mm_pCT_60W_debug/sub-007/batch_job_logs/tpc8e4bab0_f2a1_40a7_9c57_54ef6096e020.mat; delete /project/2424103.01/thalstim_simulations/thalstim_sim/data/tussim/CTX500-026-010_79.6mm_pCT_60W_debug/sub-007/batch_job_logs/tpb5b9b4be_c0ba_4372_a00e_bf852ef7b5b1.m;
} 
jkosciessa commented 1 month ago

Debugging in the k-Wave function kSpaceFirstOrder suggests that this may arise from unreasonable tissue values somewhere. A major pressure line in this instance is estimated around the neck area after only a few executed time steps (far away from the transducer plane, suggesting a computational problem). I attached some plots of different time steps, that ultimately converge in NaNs. Note that this is from the pseudoCT implementation that uses continuous tissue property values.

stepa stepb stepc

jkosciessa commented 1 month ago

This seems fixed in commit 361c869e0bee9c688d96626abfef36a77f68847b after regularizing minimum pseudoCT-based bone density to the water estimate. As also discussed here this suggests that the conversion of pseudo-HU into tissue values may not be very sensible yet and should be double-checked.

I additionally increased the time step to 10% of the limit (from 1%) in commit d83e970ae456580f4a68052584014b13eae07e9f, but the former should be the relevant change.