Surface Oversampling#
Dependencies#
Crycoat should be installed
pip install cryocat
Napari should be installed within the environment for the contour drawing:
pip install napari
The ipython environment should be installed as well for napari to run from the notebook:
pip install --user ipykernel
[ ]:
%load_ext autoreload
%autoreload 2
# Import libraries
import mrcfile
import napari
import os
import numpy as np
from cryocat import cryomotl
from cryocat import cuboid_sampling
Expected time#
All steps (apart from the contour drawing that is user-dependent) should be completed within seconds/minutes
Input data#
The data for this tutorial can be downloaded here.
In the inputs folder are all files necessary to run following codes.
In the expected_outputs folder are output files that should be produced by this analysis.
Draw contour on Napari#
In the napari GUI, a new shape layer was created, and the contour was drawn using the “add polygons” button. Drawing contours on multiple axes can enhance the accuracy of sampling points on the target surfaces.
[ ]:
mrcfolder = './inputs/'
shapefolder = './inputs/'
[ ]:
tomo_num = 40
padd_num = str(tomo_num).zfill(3)
obj_num = 1
sample_distance = 3
[ ]:
#Initialize the napari viewer:
viewer = napari.Viewer()
# Load tomogram and add it to the viewer
tomo = mrcfile.open(mrcfolder + padd_num + '.mrc')
tomo_data = tomo.data
viewer.add_image(tomo_data)
[ ]:
# opening the shape layer which is provided in inputs folder
polygons = cuboid_sampling.load_shapes(shapefolder+"040_1_shape.csv")
shapes_layer = viewer.add_shapes(polygons, shape_type='polygon')
[ ]:
# For own data one has to create the shape layer and draw the polygons
shapes_data=viewer.layers[1].data
# And save it for future work
cuboid_sampling.save_shapes(shapes_data, shapefolder + padd_num +'_'+str(obj_num)+'_shape.csv') # this file is available in the expected_outputs to continue further
# The following code will produce the same results only for the 040_1_shape.csv (provided). For newly produced shapes the outcome will be of course different.
Surface sampling from shapes#
A list with the same length as the number of shapes could serve as the input for sampling and shifting distances. To exclude both the top and bottom surfaces from oversampling, set tb_dist to 0. Use 1 to exclude the face with a higher z value, or -1 to exclude the face with a smaller z value. Record the shift_dist and tb_dist. This would used for resetting normal.
[ ]:
output_path = './outputs/'
shapes_path = './inputs/'
overSample_dist = 3
shift_dist = 4
tb_dist = 0
rm_surface = None
[ ]:
shapes = os.listdir(shapes_path)
for i,name in enumerate(shapes):
if name.endswith(".csv"):
tomo_num = int(name[0:3])
obj_num = int(name[4:5])
if type(overSample_dist) == int:
overSample_dist = overSample_dist
elif type(overSample_dist) == list:
overSample_dist = overSample_dist[i]
if type(shift_dist) == int:
shift_dist = shift_dist
elif type(shift_dist) == list:
shift_dist = shift_dist[i]
if type(tb_dist) == int:
tb_dist = tb_dist
elif type(tb_dist) == list:
tb_dist = tb_dist[i]
shapes_data=cuboid_sampling.load_shapes(shapes_path+name)
pd_points, pd_angles = cuboid_sampling.get_sampling_pandas(shapes_data, overSample_dist, shift_dist, tb_dist, rm_surface)
pd_motl = cryomotl.EmMotl()
pd_motl.fill({"coord": pd_points, "angles": pd_angles, "tomo_id": tomo_num, "object_id":obj_num})
pd_motl.write_out(output_path + str(tomo_num).zfill(3) + '_' + str(obj_num) + '_motl_sp'+ str(overSample_dist) +'.em')
Reset normal of points#
[ ]:
motl = cryomotl.Motl.load("./inputs/dist_clean20_1.em")
motl = motl.get_motl_subset(40)
motl.write_out("./inputs/dist_clean20_1.em")
[ ]:
motlfolder = './inputs/dist_clean20_1.em'
outfolder = './outputs/dist_clean20_reNorm_1.em' # available in expected_outputs
searching_sampDist = 1 # use small sample distance for searching
shifting_dist = 4 # shifting_distance for oversample points
tb_move = 0 # same as tb_dist in oversampling
binFactor = 4 # bin factor between your shape and motl_list
motl = cryomotl.Motl.load(motlfolder)
tomoObjNum = motl.df.values[:,4:6]
uniTomoObjNum = np.unique(tomoObjNum, axis = 0)
tomoNum = 0
for i in uniTomoObjNum:
lastTomoNum = tomoNum
tomoNum = int(i[0])
objNum = int(i[1])
if tomoNum == lastTomoNum:
shapesObjNum = shapesObjNum+1
else:
shapesObjNum = 1
# find the rows belongs to the shape
isRow = [i[5] == objNum for i in motl.df.values[:]]
df_motl = motl.df.loc[isRow]
# load plygons shapes
shapes_data = cuboid_sampling.load_shapes(shapefolder + str(tomoNum).zfill(3) +'_' + str(shapesObjNum) + '_shape.csv')
# The motl input of reset_normals needed to be match with the shape
df_motl = cuboid_sampling.reset_normals(searching_sampDist, shifting_dist, tb_move, shapes_data, df_motl, binFactor)
# insert the updated motl back into motl with mutiple tomograms
motl.df.update(df_motl)
cryomotl.Motl.write_out(motl, outfolder)
Calculate surface area#
You have the option to exclude the top or bottom surfaces from the calculation. If you wish to remove both the top and bottom surfaces, you can set rm_faces_all to 0.
[ ]:
pixel_size = 1.56
binning = 8
rm_faces_all = 0
shapes = os.listdir('./outputs/')
np_area = np.eye(len(shapes),3)
for n,i in enumerate(shapes):
if i.endswith(".csv"):
tomo_num = int(i[0:3])
object_num = int(i[4:5]) # not the same with motl list
padd_num = str(tomo_num).zfill(3)
shapes_data=cuboid_sampling.load_shapes(shapefolder+ padd_num +'_' + str(object_num) + '_shape.csv')
if isinstance(shapes_data, list):
# create array from the list
mask_points=np.concatenate(shapes_data, axis=0 )
else:
mask_points=shapes_data
if rm_faces_all == 0:
rm_faces = rm_faces_all
else:
rm_faces = rm_faces_all[n]
sur_area = cuboid_sampling.get_surface_area_from_hull(mask_points, rm_faces)*(pixel_size**2)*(binning**2)
np_area[n,0] = tomo_num
np_area[n,1] = object_num
np_area[n,2] = sur_area
np.savetxt('./outputs/invitro_area.csv',np_area,delimiter=',') # available in expected_outputs