TAP/TMS_forwardsolution_angle_to_midline.m at main · moritzdannhauer/TAP · GitHub

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
clear all;

TMS_COIL_ANGLE_TO_MIDLINE = 45; %degree
HairThinkness = 2; %mm

if (ispc)
    sep='\';
elseif (ismac || isunix)
    sep='/';
end

subj='ernie_headreco_no-conform';
subjects_folder='Subjects';
target_name='DLPFC';
mask='DLPFC_mask_TMS_target.nii.gz';
outputfolder='DLPFC_TMSsim';
scirunfolder='DLPFC_scirun';
SimNIBS_MRI=[subj '_T1fs_conform.nii.gz'];

simnibs_folder='/Users/rothermichk17/Applications/SimNIBS-3.2/';
root_folder='/Users/rothermichk17/Pipeline_git';
freesurfer_matlab_folder='/Applications/freesurfer/7.1.0//matlab/';

coil_name='Medtronic_MCF_B65.ccd';
addpath(char([simnibs_folder sep 'matlab']));
addpath(char([root_folder sep 'matlab']));
addpath(char(freesurfer_matlab_folder));

[brain, scalp, head_model_mesh, brain_vert_neigh]=get_meshSurfs(subjects_folder, subj, [1 2], sep);
[coord, which_pipeline] = get_target(subjects_folder,subj,mask,SimNIBS_MRI, sep, brain);

scalp_nodes=scalp.node';
[~, index, ~] = min_distance(scalp_nodes, coord);
spp=scalp_nodes(index,:);
zdir=scalp.field(:,index)';

v=[0 1 0]; %ant-post direction
ydir=cross(v, zdir);
ydir=ydir/norm(ydir);
xdir=cross(ydir, zdir);
xdir=xdir/norm(xdir);

u=zdir;
ux=u(1);
uy=u(2);
uz=u(3);
c=cosd(TMS_COIL_ANGLE_TO_MIDLINE);
s=sind(TMS_COIL_ANGLE_TO_MIDLINE);
R=[c+ux^2*(1-c) ux*uy*(1-c)-uz*s ux*uz*(1-c)+uy*s;
     uy*ux*(1-c)+uz*s c+uy^2*(1-c) uy*uz*(1-c)-ux*s;
     uz*ux*(1-c)-uy*s uz*uy*(1-c)+ux*s c+uz^2*(1-c)
     ];

xdir=R*xdir';
ydir=R*ydir';

%for visualization in scirun
Field1.node=spp';
Field1.field=xdir;
Field2.node=spp';
Field2.field=ydir;
Field3.node=spp';
Field3.field=zdir';

T=zeros(4,4);
T(1:3,1)=xdir;
T(1:3,2)=ydir;
T(1:3,3)=-zdir;
tmp=spp+HairThinkness.*zdir;
T(1:3,4)=tmp';
T(4,4)=1;

Flip_current_direction=0;

S = sim_struct('SESSION');
S.fnamehead = [subjects_folder sep subj sep 'ernie_headreco_no-conform.msh']; % head mesh
S.pathfem = [ subjects_folder sep subj sep outputfolder];
S.poslist{1} = sim_struct('TMSLIST');
S.poslist{1}.fnamecoil = 'Magstim_70mm_Fig8.ccd'; % Choose a coil in the ccd-files folder
S.poslist{1}.anisotropy_type='scalar';
S.poslist{1}.aniso_maxcond=2;
S.poslist{1}.aniso_maxratio=10;
S.poslist{1}.pos(1).distance = 0; % 4 mm distance from coil surface to head surface
S.poslist{1}.pos(1).didt = 1e6; % 4 mm distance from coil surface to head surface
S.poslist{1}.pos(1).matsimnibs=T;
run_simnibs(S);

A=T;

A(:,3)=-A(:,3); %Brainsight normal is inverted compared to SimNIBS, so flipping it here

make_brainsight_files(subjects_folder, subj, target_name, A, outputfolder, HairThinkness, Flip_current_direction, sep, 2);