In [1]:
%matplotlib inline
import nibabel as nib
import numpy as np
import matplotlib.pyplot as plt
import os
import sys
sys.path.insert(0,'../')
import surfdist as sd
from surfdist import viz, load, utils, surfdist


/usr/local/lib/python2.7/site-packages/matplotlib/font_manager.py:273: UserWarning: Matplotlib is building the font cache using fc-list. This may take a moment.
  warnings.warn('Matplotlib is building the font cache using fc-list. This may take a moment.')

In [8]:
# load surface data and other variables
cmap = 'coolwarm'
base_dir = '/Applications/freesurfer/subjects/'
surf = nib.freesurfer.read_geometry(os.path.join(base_dir, 'bert/surf/lh.pial'))
cort = np.sort(nib.freesurfer.read_label(os.path.join(base_dir, 'bert/label/lh.cortex.label')))
sulc = nib.freesurfer.read_morph_data(os.path.join(base_dir, 'bert/surf/lh.sulc'))

Calculate and display distance from central sulcus at each node:


In [ ]:
# print all label names
sd.load.get_freesurfer_label(os.path.join(base_dir, 'bert/label/lh.aparc.a2009s.annot'))

In [ ]:
# load central sulcus nodes
region = 'S_central'
src  = sd.load.load_freesurfer_label(os.path.join(base_dir, 'fsaverage4/label/lh.aparc.a2009s.annot'), region, cort)

# calculate distance
dist = sd.surfdist.dist_calc(surf, cort, src)
print np.shape(dist)

In [11]:
dist1 = dist.copy()
dist1[cort] = (dist1[cort] - np.max(dist1[cort])) * -1

In [13]:
# visualize
plot_med, ax_med = sd.viz.viz(surf[0], surf[1], dist1, bg_map=sulc, bg_on_stat=True, cmap=cmap, figsize=((15,10)))
ax_med.axes.set_axis_off()
ax_med.grid(False)
plt.show()
plot_med.savefig('fig.dist1.med.png')
plot_lat, ax_lat = sd.viz.viz(surf[0], surf[1], dist1, azim=180, bg_map=sulc, bg_on_stat=True, cmap=cmap, figsize=((15,10)))
ax_lat.axes.set_axis_off()
ax_lat.grid(False)
plt.show()
plot_lat.savefig('fig.dist1.lat.png')



In [16]:
# visualize
plot_med, ax_med = sd.viz.viz(surf[0], surf[1], dist2, bg_map=sulc, bg_on_stat=True, cmap=cmap, figsize=((15,10)))
ax_med.axes.set_axis_off()
ax_med.grid(False)
plt.show()
plot_med.savefig('fig.dist2.med.png')
plot_lat, ax_lat = sd.viz.viz(surf[0], surf[1], dist2, azim=180, bg_map=sulc, bg_on_stat=True, cmap=cmap, figsize=((15,10)))
ax_lat.axes.set_axis_off()
ax_lat.grid(False)
plt.show()
plot_lat.savefig('fig.dist2.lat.png')



In [37]:
d = np.argmin(np.vstack((dist,distX)), axis=0)
d[cort] = d[cort] + 1

In [42]:
# visualize
plot_med, ax_med = sd.viz.viz(surf[0], surf[1], d, bg_map=sulc, bg_on_stat=True, cmap=cmap, figsize=((15,10)))
ax_med.axes.set_axis_off()
ax_med.grid(False)
plt.show()
plot_med.savefig('fig.dist3.med.png')
plot_lat, ax_lat = sd.viz.viz(surf[0], surf[1], d, azim=180, bg_map=sulc, bg_on_stat=True, cmap=cmap, figsize=((15,10)))
ax_lat.axes.set_axis_off()
ax_lat.grid(False)
plt.show()
plot_lat.savefig('fig.dist3.lat.png')



In [39]:
# save colorbar
gradient = np.linspace(0, 1, 256)
for x in range(5):
    gradient = np.vstack((gradient, gradient))

fig, ax = plt.subplots(nrows=1)
ax.imshow(gradient, cmap=cmap)
ax.set_axis_off()
fig.tight_layout()
plt.show()
fig.savefig('fig.cmap.colorbar.png')

print(np.max(dist))
print(np.min(dist))


179.858349943
0.0

In [32]:
# Calculate zones on individual
# load central sulcus and calcarine
regions = ['S_calcarine', 'S_central']
src = []
for r in regions:
    src.append(sd.load.load_freesurfer_label(os.path.join(base_dir, 'bert/label/lh.aparc.a2009s.annot'), r, cort))

# calculate zones
zone = sd.surfdist.zone_calc(surf, cort, src)

In [34]:
zones = zone.copy()
print np.unique(zone)

zone[np.where(zone == 1)] = -1
zone[np.where(zone == 2)] = 1


[ 0.  1.  2.]

In [35]:
# visualize zones
plot_zones_med, ax_zones_med = sd.viz.viz(surf[0], surf[1], zone, bg_map=sulc, bg_on_stat=True, cmap=cmap, figsize=((15,10)))
ax_zones_med.axes.set_axis_off()
ax_zones_med.grid(False)
plt.show()
plot_zones_med.savefig('fig.zones.med.png')
plot_zones_lat, ax_zones_lat = sd.viz.viz(surf[0], surf[1], zone, azim=180, bg_map=sulc, bg_on_stat=True, cmap=cmap, figsize=((15,10)))
ax_zones_lat.axes.set_axis_off()
ax_zones_lat.grid(False)
plt.show()
plot_zones_lat.savefig('fig.zones.lat.png')