Skip to content

Commit

Permalink
Refactoring
Browse files Browse the repository at this point in the history
  • Loading branch information
@jcohenadad
jcohenadad committed Jan 31, 2024
1 parent 62acbeb commit 24b3f00
Showing 1 changed file with 51 additions and 59 deletions.
110 changes: 51 additions & 59 deletions data_processing.ipynb
View file
Original file line number Diff line number Diff line change
Expand Up @@ -867,66 +867,58 @@
"font_size = 14\n",
"axes=axes.flatten()\n",
"\n",
"# Plotter loop (avoiding the generation of an ungly 4D data structure)\n",
"# First, plot the anatomical image with an overlay of the mask\n",
"\n",
"# Load data\n",
"CP_anat=nib.load(f\"{subject}_acq-anatCP_TB1TFL.nii.gz\")\n",
"CP_SC=nib.load(file_mask)\n",
"CP_nTpV=nib.load(f\"{subject}_acq-CP_TB1map.nii.gz\")\n",
"\n",
"# Defining mask based on the magnitude image intensity threshold\n",
"cslice=CP_anat.shape[2] // 2 -2 #shows the SC seg best\n",
"threshold=300\n",
"mask=CP_anat.get_fdata()[xmin:xmax,ymin:ymax, cslice]\n",
"mask=np.where(mask > threshold, 1, 0)\n",
"\n",
"# Cropping anat, SC, B1+\n",
"CP_anat=CP_anat.get_fdata()[xmin:xmax,ymin:ymax, cslice]\n",
"CP_anat=CP_anat*mask\n",
"CP_SC=CP_SC.get_fdata()[xmin:xmax,ymin:ymax, cslice]\n",
"CP_nTpV=CP_nTpV.get_fdata()[xmin:xmax,ymin:ymax, cslice]\n",
"CP_nTpV=CP_nTpV*mask\n",
"\n",
"# All opacity overalys look ugly: workaround, set the anat slice to a max value where the segmentation exists\n",
"CP_anat[CP_SC>0.5]=4000;\n",
"\n",
"# Plotting anat overlayed with SC\n",
"splot=axes[0]\n",
"splot.imshow((CP_anat.T), cmap='gray', origin='lower',vmin=0,vmax=2000)#, interpolation='spline36')\n",
"splot.set_title('SC overlay', size=font_size)\n",
"splot.axis('off')\n",
"\n",
"# Then, plot each B1+ map, with an overlay of the mean and CV inside the cord\n",
"for i,shim_mode in enumerate(shim_modes):\n",
" if i==0: # Grabbing the CP mode anat to use for displaying the SC and cropping noise for the B1+ map\n",
" \n",
" # Load data\n",
" CP_anat=nib.load(f\"{subject}_acq-anat{shim_mode}_TB1TFL.nii.gz\")\n",
" CP_SC=nib.load(file_mask)\n",
" CP_nTpV=nib.load(f\"{subject}_acq-{shim_mode}_TB1map.nii.gz\")\n",
" \n",
" # Defining mask based on the magnitude image intensity threshold\n",
" cslice=CP_anat.shape[2] // 2 -2 #shows the SC seg best\n",
" threshold=300\n",
" mask=CP_anat.get_fdata()[xmin:xmax,ymin:ymax, cslice]\n",
" mask=np.where(mask > threshold, 1, 0)\n",
" \n",
" # Cropping anat, SC, B1+\n",
" CP_anat=CP_anat.get_fdata()[xmin:xmax,ymin:ymax, cslice]\n",
" CP_anat=CP_anat*mask\n",
" CP_SC=CP_SC.get_fdata()[xmin:xmax,ymin:ymax, cslice]\n",
" CP_nTpV=CP_nTpV.get_fdata()[xmin:xmax,ymin:ymax, cslice]\n",
" CP_nTpV=CP_nTpV*mask\n",
" \n",
" # All opacity overalys look ugly: workaround, set the anat slice to a max value where the segmentation exists\n",
" CP_anat[CP_SC>0.5]=4000;\n",
" \n",
" # Plotting anat overlayed with SC\n",
" splot=axes[i]\n",
" splot.imshow((CP_anat.T), cmap='gray', origin='lower',vmin=0,vmax=2000)#, interpolation='spline36')\n",
" splot.set_title('SC overlay', size=font_size)\n",
" splot.axis('off')\n",
" \n",
" # Plotting the B1+ map\n",
" splot=axes[i+1]\n",
" splot.imshow((CP_nTpV.T), cmap='viridis', origin='lower',vmin=dynmin,vmax=dynmax)#, interpolation='spline36')\n",
" splot.set_title(shim_mode, size=font_size)\n",
" splot.axis('off')\n",
"\n",
" \n",
" else:\n",
" # Load data\n",
" B1map=nib.load(f\"{subject}_acq-{shim_mode}_TB1map.nii.gz\")\n",
" B1map=B1map.get_fdata()[xmin:xmax,ymin:ymax, cslice]\n",
" B1map=B1map*mask\n",
" \n",
" # Plot\n",
" splot=axes[i+1]\n",
" im = splot.imshow((B1map.T), cmap='viridis', origin='lower',vmin=dynmin,vmax=dynmax)#, interpolation='spline36')\n",
" splot.set_title(shim_mode, size=font_size)\n",
" splot.axis('off')\n",
"\n",
" # Find the statistics for the current shim mode\n",
" shim_stats = subject_stats[subject_stats['Shim_Mode'] == shim_mode]\n",
" if not shim_stats.empty:\n",
" mean_val = shim_stats.iloc[0]['Average']\n",
" std_val = shim_stats.iloc[0]['Standard_Deviation']\n",
" cv = std_val / mean_val * 100 # Coefficient of variation in percentage\n",
" annotation_text = f\"{mean_val:.2f} nT/V\\n{cv:.2f}%\"\n",
" splot.annotate(annotation_text, (0.05, 0.95), xycoords='axes fraction', \n",
" fontsize=10, color='white', \n",
" verticalalignment='top', horizontalalignment='left')\n",
" # Load data\n",
" B1map=nib.load(f\"{subject}_acq-{shim_mode}_TB1map.nii.gz\")\n",
" B1map=B1map.get_fdata()[xmin:xmax,ymin:ymax, cslice]\n",
" B1map=B1map*mask\n",
"\n",
" # Plot\n",
" splot=axes[i+1]\n",
" im = splot.imshow((B1map.T), cmap='viridis', origin='lower',vmin=dynmin,vmax=dynmax)#, interpolation='spline36')\n",
" splot.set_title(shim_mode, size=font_size)\n",
" splot.axis('off')\n",
"\n",
" # Find the statistics for the current shim mode\n",
" shim_stats = subject_stats[subject_stats['Shim_Mode'] == shim_mode]\n",
" if not shim_stats.empty:\n",
" mean_val = shim_stats.iloc[0]['Average']\n",
" std_val = shim_stats.iloc[0]['Standard_Deviation']\n",
" cv = std_val / mean_val * 100 # Coefficient of variation in percentage\n",
" annotation_text = f\"{mean_val:.2f} nT/V\\n{cv:.2f}%\"\n",
" splot.annotate(annotation_text, (0.05, 0.95), xycoords='axes fraction', \n",
" fontsize=10, color='white', \n",
" verticalalignment='top', horizontalalignment='left')\n",
"\n",
"plt.tight_layout()\n",
"plt.subplots_adjust(wspace=0.1, hspace=0, right=0.9)\n",
Expand Down

0 comments on commit 24b3f00

Please sign in to comment.