diff --git a/data_processing.ipynb b/data_processing.ipynb index 2385dba..d176545 100644 --- a/data_processing.ipynb +++ b/data_processing.ipynb @@ -439,6 +439,7 @@ "\n", "# Adjust the layout so labels and titles do not overlap\n", "plt.tight_layout()\n", + "plt.savefig(os.path.join(path_results, 'fig_gre_csf-sc_ratio.png'), dpi=300, format='png')\n", "plt.show()" ] }, @@ -734,6 +735,7 @@ "\n", "# Adjust the layout so labels and titles do not overlap\n", "plt.tight_layout()\n", + "plt.savefig(os.path.join(path_results, 'fig_b1plus.png'), dpi=300, format='png')\n", "plt.show()" ] }, @@ -843,6 +845,14 @@ "\n", "# Select subject to show\n", "subject = 'sub-05'\n", + "os.chdir(os.path.join(path_data, subject, \"fmap\"))\n", + "file_mask = f\"{subject}_acq-anatCP_TB1TFL_mask-shimming.nii.gz\"\n", + "\n", + "# Load the statistics from the CSV file\n", + "stats_df = pd.read_csv(os.path.join(path_results, 'stats_b1plus.csv'))\n", + "\n", + "# Filter for the specific subject\n", + "subject_stats = stats_df[stats_df['Subject'] == subject]\n", "\n", "# Defining crop limits for resulting figure\n", "xmin = 20\n", @@ -859,55 +869,58 @@ "font_size = 14\n", "axes=axes.flatten()\n", "\n", - "# Plotter loop (avoiding the generation of an ungly 4D data structure)\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", + "# First, plot the anatomical image with an overlay of the mask\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", + "# 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('Anat', 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", + " # 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", @@ -923,6 +936,7 @@ "# cbar_ax = fig.add_axes([0.95, 0.5, 0.04, 0.4])\n", "# cbar = plt.colorbar(im, cax=cbar_ax)\n", "cbar_ax.set_title('nT/V', size=12)\n", + "plt.savefig(os.path.join(path_results, 'fig_b1plus_map.png'), dpi=300, format='png')\n", "plt.show\n" ] },