From 422eec3a8c0eefe9a790c9fc7beec9b8381fdb83 Mon Sep 17 00:00:00 2001 From: cgohil8 <53237863+cgohil8@users.noreply.github.com> Date: Fri, 1 Dec 2023 13:27:57 +0000 Subject: [PATCH] RHINO example without visualisations + AAL parcellation. (#240) --- examples/oxford/no_viz_pipeline/README.md | 4 + examples/oxford/no_viz_pipeline/do_rhino.py | 100 ++++++++++++++++++ examples/oxford/ssp_pipeline/4_coregister.py | 2 +- .../aal_cortical_merged_8mm_stacked.nii.gz | Bin 0 -> 36166 bytes 4 files changed, 105 insertions(+), 1 deletion(-) create mode 100644 examples/oxford/no_viz_pipeline/README.md create mode 100644 examples/oxford/no_viz_pipeline/do_rhino.py create mode 100644 osl/source_recon/parcellation/files/aal_cortical_merged_8mm_stacked.nii.gz diff --git a/examples/oxford/no_viz_pipeline/README.md b/examples/oxford/no_viz_pipeline/README.md new file mode 100644 index 00000000..b8f9bd9d --- /dev/null +++ b/examples/oxford/no_viz_pipeline/README.md @@ -0,0 +1,4 @@ +Pipeline without any visualisations +----------------------------------- + +Sometimes there are errors when generating plots. In this pipeline we use OSL without any visualisations. diff --git a/examples/oxford/no_viz_pipeline/do_rhino.py b/examples/oxford/no_viz_pipeline/do_rhino.py new file mode 100644 index 00000000..f08d817b --- /dev/null +++ b/examples/oxford/no_viz_pipeline/do_rhino.py @@ -0,0 +1,100 @@ +"""RHINO without any visualisations. + +""" + +# Authors: Chetan Gohil + +import numpy as np + +from osl.source_recon import setup_fsl, rhino + +setup_fsl("/opt/ohba/fsl/6.0.5") # this is where FSL is installed on hbaws + +# Directories +preproc_dir = "output/preproc" +anat_dir = "path/to/smri/dir" +coreg_dir = "output/coreg" + +# Files ({subject} will be replaced by the name for the subject) +preproc_file = preproc_dir + "/{subject}_tsss_preproc_raw.fif" +smri_file = anat_dir + "/{subject}/anat/{subject}_T1w.nii" + +# Subjects to coregister +subjects = ["sub-001", "sub-002"] + +def fix_headshape_points(src_dir, subject, preproc_file, smri_file): + filenames = rhino.get_coreg_filenames(src_dir, subject) + + # Load saved headshape and nasion files + hs = np.loadtxt(filenames["polhemus_headshape_file"]) + nas = np.loadtxt(filenames["polhemus_nasion_file"]) + lpa = np.loadtxt(filenames["polhemus_lpa_file"]) + rpa = np.loadtxt(filenames["polhemus_rpa_file"]) + + # Remove headshape points on the nose + remove = np.logical_and(hs[1] > max(lpa[1], rpa[1]), hs[2] < nas[2]) + hs = hs[:, ~remove] + + # Remove headshape points on the neck + remove = hs[2] < min(lpa[2], rpa[2]) - 4 + hs = hs[:, ~remove] + + # Remove headshape points far from the head in any direction + remove = np.logical_or( + hs[0] < lpa[0] - 5, + np.logical_or( + hs[0] > rpa[0] + 5, + hs[1] > nas[1] + 5, + ), + ) + hs = hs[:, ~remove] + + # Overwrite headshape file + np.savetxt(filenames["polhemus_headshape_file"], hs) + + +# Main RHINO code + +for subject in subjects: + + # Input files + preproc_file_ = preproc_file.format(subject=subject) + smri_file_ = smri_file.format(subject=subject) + + # Extract fiducials and headshape points + filenames = rhino.get_coreg_filenames(src_dir, subject) + rhino.extract_polhemus_from_info( + fif_file=preproc_file_, + headshape_outfile=filenames["polhemus_headshape_file"], + nasion_outfile=filenames["polhemus_nasion_file"], + rpa_outfile=filenames["polhemus_rpa_file"], + lpa_outfile=filenames["polhemus_lpa_file"], + ) + + # Clean up headshape points + fix_headshape_points(src_dir, subject, preproc_file_, smri_file_) + + # Extract surfaces from the structural + rhino.compute_surfaces( + smri_file=smri_file_, + subjects_dir=src_dir, + subject=subject, + include_nose=False, + ) + + # Coregistration + rhino.coreg( + fif_file=preproc_file_, + subjects_dir=src_dir, + subject=subject, + use_headshape=True, + use_nose=False, + ) + + # Forward model + rhino.forward_model( + subjects_dir=src_dir, + subject=subject, + model="Single Layer" + gridstep=8, + ) diff --git a/examples/oxford/ssp_pipeline/4_coregister.py b/examples/oxford/ssp_pipeline/4_coregister.py index d6e0caa0..02a11c34 100644 --- a/examples/oxford/ssp_pipeline/4_coregister.py +++ b/examples/oxford/ssp_pipeline/4_coregister.py @@ -17,7 +17,7 @@ from osl import source_recon, utils # Directories -preproc_dir = "/output/preproc_ssp" +preproc_dir = "output/preproc_ssp" anat_dir = "path/to/smri/dir" coreg_dir = "output/coreg" fsl_dir = "/opt/ohba/fsl/6.0.5" # this is where FSL is installed on hbaws diff --git a/osl/source_recon/parcellation/files/aal_cortical_merged_8mm_stacked.nii.gz b/osl/source_recon/parcellation/files/aal_cortical_merged_8mm_stacked.nii.gz new file mode 100644 index 0000000000000000000000000000000000000000..4bce3d5663df716c61e232326c49681bbfc9a04a GIT binary patch literal 36166 zcmeHwdtB1@{(tjf=U6*CYire6TjkcVwmEa@%*^Zg?rg41ee5PFB5kJR5Ghef5FBUU z?{>B>HcO2}==(Ku9ThSJ@`8AFwoD}`O%YIVOC>_2L?J=o_oj2s{Jh;c+gZO|s`-x( z^>cl{F3;Ef{aLf$%S-S90qx;>`|V6@Dlt1JBNhBYNXy2h?ZR#*5U`}2w^Q+HyEbHH zWGomP{OiHz6O$GW1$Wo(So{0To&U2isbog4Ugesxiy{l{YGR?wXliv#sX27JWX&2| zmq;>W`D`s~6eYXb z!A?I8{r4x7m%Ei@KdI9E;6xY-@+|GjUDF}NRH4q+VQmf8nE0)Tr!dGXc85#AF=Kd@ z;#gg{nxcORargm`$i>42(^so?Sa4VuY_I==&2j6rFN`8Ll00M*-c*F~gL!%RVL*0# z>5dK7LGmONSKDtGEb5WKbQ$ZDlny3eB;Cdt zE=>&N#Bw&`xby^@q)AXJfs^#RSLXQ};IHKnweMSqeRU z&#AWT5Ym(yUyYiuROn{biHmg&XP%qmm=GC#q&7)gCJ z7otx%3N$KzaWP5r!L#t=PP}o#F{9S-%ysV~)r@g)CQ^>5^Yj_B@1^#fAZKK&8>Zu{QQ6ML>UF13vwxA+2`w?aVb-9myf~#J z|A;Q`T~jbwbZS=6cC0@vDL-O3HD-O2;$uZ;3=U6jik9s-s5Lw~myhnDZy5Y5{i{!fnl zgDyKC+Y_NZ)tMeN6!)TcG!RN~i z;to3}XZ;*meIdqMH!)n>I0`?AvCbWZ#16KSGmO2XCm)`3)!pG9+;=*-^RXQp6!FDZ zP4cEV&PMt5K(1PaAJZTI8l2p}G%aM06ev}E0jf3E_VQQt2h^$YjpfKwc>b3tDcg9% z2FapOWWce2L0-N)(j~M~b*4jNNF9-v&%3Jso*Fwr_}3U0BV&b=5~va;HE2I{t|tD3 zV?eFq2^4nP`^IC?u|rFChyS`-aV+}S0RKFQ;;mtZ)i6-Lu(*Nk+v6Kz$!d&~?26&l zM^h#BJUsHGmhvbWp+o2pS>7deL#dj1u0hm>Ov=ADBe+EX6gyMdz3ojf@rG6{JWrb5 zM-z7R<60t`Ia6?wyKOt;=0LfDiLYh$cBjUY<#V7D^~iZTT5OIIo zKFungNXovMX*&>-YAJm@m6J?lhM1L;TmKCzazL7Ee>Yw#x9mBOOXc0P9f(__c+x&> zyBbcPPjXOju)!6K=+HelRB+_XvLKJ5xR$OQ#OG_(? z#ODpAZ2e`8e>fcI%3_j>r3U_Y!Yad;@bEh>>a-^VE?0^vl=4Qo@E#E{ghu7M}XroMfF&&+0GB7oeyg+gYDW-qvv^28W)}S ziTli+_aUmXU)mXl8f~hg*$$k>?tc${!hAGv$BKI9L$^PZF4l&#w{kXG zci_9Elj?xK@gY{efc4f1&>@9iK?s7k`gbS?JKqsH7xn3fdkC#oUZ70i+YTEXiN5~4Fq9m&( z7kH^uj@I5J!V zoRaH`r%yoi`M%G*&q(+qFUJUMo9YE+U$7YXin{*w3Y@jmMxb@aJc`V_wviqjZV4Gx zHve7N`_})+Z#9qJkYD~CDloK}F_Mk(%^&E*@{UHN(j8ZBN@68fQNpK@IOA3O$nNlV z(3*$0oNEuP6Q5N~!ksD$7;y}iw%E_FN0%U+%G6iFY8Y15R!1U6L{fwG8%LGmvvfRf z@(7Yqo0E2itsk*nTYZ(H-|O3MFKW#)F8i<^8$s-t#KU(}FPtFLv+bXUEbZ1M$P&{W zN73SmXl8Do>6_n|k!^gsl@YklgLBW>Aw4Z-H4R-O7+ z>7G+jhUs8&aiBF5Hoh&KDEst8qD^2ekPwBdkhO%NkcBE)8E9xK&!S2;`6Jug_9c7k z4z+W_EFW3e`m(SDowGN|UWJ-)u4+7jBs0(FFzL0k1hT#)Q(p#vUliwk$iTKdhNLNd z!(meNgf5;%ipEPOx3Mcp)ZWQ$npX3Gw!s>TYQH|qPxebE93~r1uTTYo);MkFzs^$tn_xW9U4ih(Jajd9Tgjwb_+LEi=;9 zc{+%lOEzxojQIVz6bd$q4gU(%AR zP2_-9$J*7AXpvYex7ouXp=S9WCZ$abY${9+7Z(?{R9-z|p|9GO1mRtqpIaY9s=o<9Aa~X1ryMO#+3M(9)L^SqAFM_S%Z+nUvX?vB*N@r% z2pV)S#+!P5NkPG9t??&O+1lRikT1dzc!Y?26f&iuadO+ZK&EQcmn4U_lOUen(0GKP z)}ez`c0rRK+4To48d@C^x1f6=Ow`YRj^+98)oxNQw?Gp?6gFiOChfl-UPf zvjXGzX|4%VEkM~sbB8V}Jyr6pHiCU!gH2%hIFBlsJPHfoPSC^{ve{7_XvzH5y)MA_ z97n&E;PR)=BsWn{Vb3<#N zj!RYqT2^0|Sf@9ro<%%d?l&|eG&*joM`7dlq!0LitSxaSi7K?HEOj?t{Q4q>^|9S9 zyAybKmP})bOIU>s;D~nNdCS*?&$6MweNKp@cujNg%z(;wlf0WJz)^gzK8NrK_C zJVTGaavcjG`0Ay9Mo-jNTP4i`g)>y-@%8f_@Q8k}ffG3yYJ1}CxxHgQn_pyqaJaBN ztkLzcA>mcYaMQDFCJY5piD8Q0t7vJ)KQFiSVsLarQeVa6RUV&C+ zud4n!^R(|_K%96$oQZ^P+;%{n-#M}!{YM?9INq?+HwhY~sQ_-((+(kxsB#F6vK~le z=7D!sD8Jh>6F|YN?y?{r&x`Y`LteP2HA=R(!~6v1iIehUeva<7j_kowqeIn#n54OZ zY13XdwElW;47WVMb$aq~l5?ASozSWlhQ`l^K*!dCwqf1&@<<;?u^H9BX`{GUed=Q2 z7V$s+IzO-}SJRWqOO~VGR-Y6`sPnEO+)^+dQ>Iba9HnGr}Lg(m)m1z2_|%g zQ}Gk#=7McM8)VJH!UBd4%sdWO08~{8#n_s1s0XMYeZtYwry-iGsZbss2?O(J`W@8R zDQwGThVLD3vCv6m{7Qro_hh@fAN$K*&XSx7me4o`!!~D z(&Ggd6|sHYtP%9yQwkC}kt*O^;*BLHDu_-pW9fCoSNvz=%c2JxTr{5cs_DwodFc8+ z5YSpf$4ZjX-aC*t%DpU2-DS2grLd?kg%D?hUa(E^_tY3&y#Hz;m?XdYL_zC2wgc1CgF4z@v`^YM z%b1%Ax7mDjyc@@t=8H$Eg0Ht!U|1i$pq)w;K63yx_VKBQXyB(k&C!8 z6G$&hNo8LisR`m#X;IY8X7bHY{`-P?cKr}aawn=a3Pv{a`P-pX40xB-OD{BrP51Yh zO>sr;_DK6W;nd@35%4^MY>kL2(R)$*qhzlUXez1zw8HM{qHYJMtB24qYKX*O>54hp z+?i_7R&+C#p1gPLh+&Osi}Shg^I@$WG~ZJ?8o3#s`;k3*VcIW0p=j?Q$O*KVY7L+~ zD8Ej^5B8;qY&I6swwpkP?@$bHmzuhbcal!Tdr= zdNqjy>LohvSz3~fz}bb&b$WWW@ZU3!gi2^ss6chcI_2a!C^4_2mH+n)h4PAv`qI*@ zflqgQIkD;N1LD$8Zwg|()8ZbQ`SR<$0g8~ zyKNuWKJc>*prSIO?@G1`)scg$lHlGIi=e&)1sLILukz)`@&&$ewx{2QJE2xG zSpI>ET0uKs-#L)_FdjO(fUz7rvK0GYaO=#xTF-`q--=`fE6! zedZF~8H-J?*-Hf(pAFCr&wT)Nn@Q-^cV-%zw>>PWnR$H+LKJT4frl;1NfSEbY)OTW zL)*A`pZO4!&-N9Oh4*70M*+9Z(Knq1=u)#C1>$>n6!q4lpzhAU=!0PrXedKQZlpCyRIG@-1zOu~FfbU*uM5U7Z7R=6 zB;$wNjXa&#s{8{vA}N5>oT{}L=JdcEzk}yBoXni=c8+ZoBQkVf9JsY0#x6}!?PiI@k9Z<@No#RJ4bdf(($Dz-Lh3oe=)*h8I^ zI%CY&*p1rPQ4t&~C%1wvpq5R6IJr$P|Smp?e{vNPq8`O*pkMxS`R{nO`7XRnmV z%Q#`99|YdgK@Q*S)nSs?6rEN5r|Cm)}53t-#@tvs9&YPmv%q8E=bUUqP%bN z;Ya+BY9X-6pGo=Dp86N%<9yhnF6n=LC%aITOud-I$EWZTsrzv_d5h+qhw;$Mx8;GJ zM1trmp>LHxN!=?#*9rVFOttu?&{tq3*o6se;rQf7Tn?n}KOa?plGnKP7{hKf#0ID} zHG-OpTj6j&Y&jUAm9>*O)(o4>cV7#D$~D?tXSI>0d+=EZ>ll(nZd)dz7R_u3C#_cA zkfzkp>M-(SP4LDnHsF{s`MwIjKk1qLDq`MA5`BEmmOkJNqP}zNb-PIdEa#LSIZig~ zbR$4wUK0KTvJ7nj))9qIN9Q*MpzW%r^!OU>SPsG~2KT4{4ZrEaw@%?@r7~rSp|n-` zNvm{P3~jozWHZ<%lvUa)B^@m#$0}G|U;rStX_Onok=y6&Y9?QPu5)^01(P62!lTNM zz)=+>J--}gg5e>%vBsBQ7l2{b(Fhm{q~t?M9wF9rVQ)1!As$x6?gL$a`X>l=u+(pD?d9D91@EE;|?7#A}#c3goYpVE!})(*K_Z8%UZ zU@9l_+}a&KC%w}TMo*q(wT5d3kN*-?R`}>StuQHrlbdAkE@_mj(Iv3jy2l5Iqm4m@ z<7}=p(sW}Hq!p~E_LPv4wgk%0P#SQQ(SjV}-*EcO-IrmbC;^?AT{Ba$3aq|s>@W7W zdNkIWD_A78cdQ2qqODq%K|`@uW915*D+{=7IBS|~(hs8}K^a6( zudqjt={lpWtm#BsfIKt5feh@+;MG9nogg1!_5~r@+!yDtm!8wIo8@rYE$Y2C3)sO| z-sp%k)g9wn=gfV1ExCm`F~wqXtETPfCkq9@GcNK!yA~o&!JndJB}_llQFh3b!0Q?k zciVy%vI4R8!BdJBIzeDGWf_#thw6|UO){*&8pi;8Z4!@I;;V2nBbh?{qK|9UoPPt7 zh5|+=!przh=6YRE$*7C%>-q;n8B;VU=zAt)+PSrg-cY{3q zgAElQ_5|t#CB2X(UBGDfK!vzYnQx+R9fff(8SSmpNb!<`6o_z}!S1D*>%U@& zJA^B>x&8TwRY=VC-t!PJ(z|Q{>Nk6WaOi3Xs+U|_4oq{h=*?=Jtgp3{o72z)l9?5a zn~=H5GlD4DDkhpI$rI;Q%hoY+34qaEi;&D&O+Sl}=6^c-cF!JKd;pEj8! z93t}7V7CZb_T6)cVFtAfvxT0HqsgD~)BcU99A$#37Iyxxtq82U1U`{6{wg-TJRPX4 zc-<0_RKboI&!V*_sy~8(0ddjYvq(sS@6k8hnK(p8gDupbpc&{3=^UeI+!%%kouewj zxD3p&_#UMKWJvB(mAVsY8|Z^XeCUs`;gwt0!FGkx2QjE{#}WVoMECmV+WuVh{7{_HN|-P z>Uf-%)xK=F&OH)qToI^)`g$NWhOE0(#fm(E$q#C5%v>=ZH z{77qp(sVi}ZCunt^{w_bi{G+PpTjWi#)hmBOY0oouJp_vzmU1xra8=&nzTVP6Sp4-!nmjY5vm2FOH+RpKzsq3e7upz_iwnmf}zPlDY(SbG9Eza#!b# zG(Sn@ykGSrdHNI(ny^bvKER&Z!*Mg=)X}NOV~Z>-Od~zV^f(8k0t1uA5a7B1d~t|X zv@)qR7^DmAeYo=JFU=XL&60hsFM$BpdMS;p= zFnN{$8r&Y%K}Cn=#F_A_4YNz=f63cbBpXOABKK{VvrM&ukwE@>1$w+by>3$4n10cZQ-NvJWETbD zt8f|I+nqa#<@lsov2KfhiCQ>)~^TFbp9n9h82p4k4W@z#FU$tC2FA0={ja@6~}F znlVhO$vS4!GeO?o+DgnDp_phKOk_AI!ro#>7@ccNcRO4Tz>xhr0bW_zQBP{ADnwl` z;tg+CyeyO2m2*+4QEFDmOJ$SN)4Z3XkQojt*vR0+2@T(0G3^`!4GQw{<92y#vfM6f zij#w4Uq7N#{_mR`eVC5r7UdIg z7U@bo;K2b693$PC`DUe}yJ{9AZ!>D_sY~?VW|`DigX`fo)R|-k>-|FRd%<%YU|9FR z;|1qW{uXU$6S%umol>D_D@e*?DFNuFoLLhHgs3gZGn>sakJ*$*4E*L%r@oV=g}*$` z_${zppcEO}Kz|E_H;Q0Pp}?A(aT?A)K0A~gGUv;l3&0n$S@NP_5!{mjJ^ou>eg4pH z&=u4AF=){5KGG^J$b!nAp!eZ6BwWbfT`8pb+Hb7(ke zTvS;CMWRDTK&V<)hnW^#(v-y^8p#;)r5Ex~PF+m0b=SX$p6OO}S3=;61F5K@?qFSF z3P}r6{`nlZpYTKQt?040r6aSu(u$54O_LtJQ~G4PyVN`Cyvr zVvaBu*%CfP0&=qo#K4);%0Qlcyvqr)8Pl3MlUDZ#?LpM&T2S)G_v44eyBGXFj3Q*k literal 0 HcmV?d00001