Skip to content

Commit

Permalink
removed uses of arraydata() (#93) + some cleaning
Browse files Browse the repository at this point in the history
  • Loading branch information
denisri committed Jan 12, 2024
1 parent cb557a8 commit da9f0dc
Show file tree
Hide file tree
Showing 7 changed files with 63 additions and 79 deletions.
10 changes: 5 additions & 5 deletions pyaims/bin/AimsCiftiSubMatrix.py
Original file line number Diff line number Diff line change
Expand Up @@ -83,19 +83,19 @@ def main():
if gyrus in labels_inv_map[i]]
if len(gyrdef) != 0:
tex, gyrus_num = gyrdef[0]
gyrus_indices = np.where(gyri[tex][0].arraydata() == gyrus_num)[0]
gyrus_indices = np.where(gyri[tex][0].np == gyrus_num)[0]
gyri_vert_l.append((tex, gyrus_indices))
else:
try:
gyrus_num = int(gyrus)
except ValueError:
print("Gyrus %s not found" % gyrus)
continue
texs = [gyrus_num in roi[0].arraydata() for roi in gyri]
texs = [gyrus_num in roi[0].np for roi in gyri]
for tex in len(texs):
if texs[tex]:
gyrus_indices = np.where(
gyri[tex][0].arraydata() == gyrus_num)[0]
gyri[tex][0].np == gyrus_num)[0]
gyri_vert_l.append((tex, gyrus_indices))

# concatenate by texture/surface
Expand All @@ -115,7 +115,7 @@ def main():
print('row_indices:', len(row_indices))
total_vertices = sum([len(troi[i][0]) for i in range(len(troi))])
print('total_vertices:', total_vertices)
col_vertices = [(i, np.where(troi[i][0].arraydata() != 0)[0])
col_vertices = [(i, np.where(troi[i][0].np != 0)[0])
for i in range(len(troi))]
col_indices = np.hstack([c.getIndicesForSurfaceIndices(1, i, vert)
for i, vert in col_vertices])
Expand All @@ -134,7 +134,7 @@ def main():
print('row:', row, ':', i, '/', len(row_indices), 'dense:',
gyr_mat.isDense())
mat.readRow(row)
row_data = mat.getRow(row).arraydata()
row_data = mat.getRow(row).np
#row_data[np.isnan(row_data)] = 0
row_data[np.isinf(row_data)] = 0
mat.freeRow(row)
Expand Down
2 changes: 1 addition & 1 deletion pyaims/python/soma/aims/meshSplit.py
Original file line number Diff line number Diff line change
Expand Up @@ -93,7 +93,7 @@ def meshSplit2(mesh, tex, graph, voxel_size=None, tex_time_step=None):
else:
tex_tstep = tex
# labels in the texture
labels = (numpy.unique(tex_tstep[0].arraydata())).tolist()
labels = (numpy.unique(tex_tstep[0].np)).tolist()
#if 0 in labels:
#labels.remove(0)

Expand Down
43 changes: 21 additions & 22 deletions pyaims/python/soma/aims/texturetools.py
Original file line number Diff line number Diff line change
@@ -1,7 +1,6 @@
#!/usr/bin/env python
# -*- coding: utf-8 -*-


from soma import aims, aimsalgo
from soma.utils.csv_utils import dict_to_table
import numpy
Expand All @@ -17,9 +16,9 @@ def mergeLabelsFromTexture(tex, labels_list, new_label):
otex: labeled texture with merged regions
"""
otex = aims.TimeTexture_S16()
tex_ar = tex[0].arraydata()
tex_ar = tex[0].np
otex[0].assign(tex_ar)
otex_ar = otex[0].arraydata()
otex_ar = otex[0].np
for i in labels_list:
otex_ar[otex_ar == int(i)] = new_label
return otex
Expand All @@ -37,8 +36,8 @@ def extractLabelsFromTexture(tex, labels_list, new_label):
otex[0].reserve(tex[0].nItem())
for i in range(tex[0].nItem()):
otex[0].push_back(0)
tex_ar = tex[0].arraydata()
otex_ar = otex[0].arraydata()
tex_ar = tex[0].np
otex_ar = otex[0].np
for i in labels_list:
otex_ar[tex_ar == int(i)] = new_label
return otex
Expand Down Expand Up @@ -69,7 +68,7 @@ def connectedComponents(mesh, tex, areas_mode=0):
with area = 16.5 and 6.0 respectively, areas are in square mm
"""
# create a numpy array from aims object
dtex = tex[0].arraydata()
dtex = tex[0].np

# number of vertices
nbvert = len(mesh.vertex())
Expand Down Expand Up @@ -97,7 +96,7 @@ def connectedComponents(mesh, tex, areas_mode=0):
otex[0].assign((dtex == label))
label_cc = aimsalgo.AimsMeshLabelConnectedComponent(mesh, otex, 0, 0)
# transform aims.TimeTexture_S16 to numpy array
label_cc_np = label_cc[0].arraydata()
label_cc_np = label_cc[0].np
step_cc[label-1].assign(label_cc_np)

if areas_mode:
Expand Down Expand Up @@ -130,19 +129,19 @@ def remove_non_principal_connected_components(mesh, tex, trash_label):
-------
out_tex: label texture
"""
t0 = tex[0].arraydata()
t0 = tex[0].np
t0 += 1 # 0 is a real label
conn_comp, areas = connectedComponents(mesh, tex, areas_mode=True)
t0 -= 1
dtype = tex[0].arraydata().dtype
dtype = tex[0].np.dtype
out_tex = aims.TimeTexture(dtype=dtype)
out_tex[0].assign(numpy.zeros(tex[0].size(), dtype=dtype))
out_arr = out_tex[0].arraydata()
out_arr = out_tex[0].np
out_arr[:] = trash_label
for label in conn_comp.keys():
comps = conn_comp[label]
largest = numpy.argmax(areas[label + 1]) + 1
comp_arr = comps.arraydata()
comp_arr = comps.np
out_arr[comp_arr==largest] = label
return out_tex

Expand All @@ -156,8 +155,8 @@ def meshDiceIndex(mesh, texture1, texture2, timestep1=0,
return
"""
tex1 = texture1[timestep1].arraydata()
tex2 = texture2[timestep2].arraydata()
tex1 = texture1[timestep1].np
tex2 = texture2[timestep2].np
if labels_table1 is not None:
tex1 = numpy.array([labels_table1[x] for x in tex1])
if labels_table2 is not None:
Expand Down Expand Up @@ -204,7 +203,7 @@ def average_texture(output, inputs):
for fname in inputs:
tex.append(aims.read(fname))
# make a 2D array from a series of textures
ar = numpy.vstack([t[0].arraydata() for t in tex])
ar = numpy.vstack([t[0].np for t in tex])
# replace the negative values by positive integers
if len(ar[ar == -1]) != 0:
tmp_label = numpy.max(ar) + 1
Expand Down Expand Up @@ -284,17 +283,17 @@ def vertex_texture_to_polygon_texture(mesh, tex, allow_cut=False):
* new mesh with possibly split triangles
It only works for meshes of triangles.
"""
dtype = tex[list(tex.keys())[0]].arraydata().dtype
dtype = tex[list(tex.keys())[0]].np.dtype
poly_tex = aims.TimeTexture(dtype=dtype)

if allow_cut:
out_mesh = mesh.__class__(mesh)

for t, tex0 in tex.items():
tdata = tex0.arraydata()
tdata = tex0.np
ptex0 = poly_tex[t]
ptex0.resize(len(mesh.polygon(t)))
poly_labels = ptex0.arraydata()
poly_labels = ptex0.np
if allow_cut:
added_vert = {}
vertex = out_mesh.vertex(t)
Expand Down Expand Up @@ -395,7 +394,7 @@ def mesh_to_polygon_textured_mesh(mesh, poly_tex):
"""
out_mesh = mesh.__class__()
out_tex = poly_tex.__class__()
dtype = poly_tex[list(poly_tex.keys())[0]].arraydata().dtype
dtype = poly_tex[list(poly_tex.keys())[0]].np.dtype
for t, tex0 in poly_tex.items():
print('t:', t)
overt = out_mesh.vertex(t)
Expand All @@ -406,8 +405,8 @@ def mesh_to_polygon_textured_mesh(mesh, poly_tex):
opoly.assign(mesh.polygon())
otex = out_tex[t]
otex.assign(numpy.zeros(mesh.vertex().size(), dtype=dtype) - 1)
#otex_arr = otex.arraydata()
tex_arr = tex0.arraydata()
#otex_arr = otex.np
tex_arr = tex0.np
added = {}
for p in range(len(mesh.polygon())):
plabel = tex_arr[p]
Expand All @@ -427,7 +426,7 @@ def mesh_to_polygon_textured_mesh(mesh, poly_tex):
poly[i] = vb
otex.data().append(plabel)
added[(v, plabel)] = vb
#otex_arr = otex.arraydata()
#otex_arr = otex.np

out_mesh.updateNormals()
return out_mesh, out_tex
Expand Down Expand Up @@ -529,7 +528,7 @@ def clean_gyri_texture(mesh, gyri_tex):
if areas_measures[label].size != 1:
wrong_labels.append(label)
for label in wrong_labels:
cc_tex_label = cc_tex[label - 1].arraydata()
cc_tex_label = cc_tex[label - 1].np
areas_measures_cc = areas_measures[label]
cc_nb = areas_measures_cc.size
for l in range(1, cc_nb):
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -9,10 +9,7 @@
# and INRIA at the following URL "http://www.cecill.info".
############################################################################

from __future__ import print_function

# system import
from __future__ import absolute_import
import numpy

# Soma import
Expand Down Expand Up @@ -41,8 +38,8 @@ def draw_sphere(mesh, longitude, latitude):
a spherical triangulation of the subject of its cortical hemisphere,
projected on a sphere
"""
lon = longitude[0].arraydata()
lat = latitude[0].arraydata()
lon = longitude[0].np
lat = latitude[0].np

nmesh = numpy.asarray(mesh.vertex())

Expand Down Expand Up @@ -162,8 +159,8 @@ def resample_mesh_to_sphere(

# ...
# multiply latitudes by 2 to get same amplitude on both coords
latitude = aims.TimeTexture(latitude[0].arraydata() * 2)
slat_tex = aims.TimeTexture(slat_tex[0].arraydata() * 2)
latitude = aims.TimeTexture(latitude[0].np * 2)
slat_tex = aims.TimeTexture(slat_tex[0].np * 2)
interpoler = aims.CoordinatesFieldMeshInterpoler(
mesh, sphere, latitude, longitude, slat_tex, slon_tex)

Expand Down Expand Up @@ -223,8 +220,8 @@ def resample_texture_to_sphere(

# ...
# multiply latitudes by 2 to get same amplitude on both coords
latitude = aims.TimeTexture(latitude[0].arraydata() * 2)
slat_tex = aims.TimeTexture(slat_tex[0].arraydata() * 2)
latitude = aims.TimeTexture(latitude[0].np * 2)
slat_tex = aims.TimeTexture(slat_tex[0].np * 2)
interpoler = aims.CoordinatesFieldMeshInterpoler(
mesh, sphere, latitude, longitude, slat_tex, slon_tex)

Expand Down Expand Up @@ -370,8 +367,8 @@ def refine_sphere_mesh(init_sphere, avg_dist_texture, current_sphere,
- init_sphere.vertex()[
init_sphere.polygon()[0][0]]).norm()

init_lat = aims.TimeTexture(init_lat[0].arraydata() * 2)
current_lat = aims.TimeTexture(current_lat[0].arraydata() * 2)
init_lat = aims.TimeTexture(init_lat[0].np * 2)
current_lat = aims.TimeTexture(current_lat[0].np * 2)
interpoler = aims.CoordinatesFieldMeshInterpoler(init_sphere,
current_sphere,
init_lat,
Expand Down
17 changes: 7 additions & 10 deletions pyaimsalgo/python/soma/aimsalgo/meshwatershedtools.py
Original file line number Diff line number Diff line change
Expand Up @@ -19,20 +19,17 @@
"""


#----------------------------Imports-------------------------------------------
# ----------------------------Imports------------------------------------------


# import system
from __future__ import absolute_import
from __future__ import print_function
import numpy
import sys

# soma import
from soma import aims, aimsalgo
from six.moves import range

#----------------------------Functions-----------------------------------------
# ---------------------------Functions-----------------------------------------


def _clean_after_merge(labels, parents, valid):
Expand Down Expand Up @@ -160,10 +157,10 @@ def watershedWithBasinsMerging(
# labelWt: array of shape (number of vertices)
# labelling of the vertices according to their basin

data = tex[0].arraydata()
idxW = idxWt[0].data().arraydata()
majorW = majorWt[0].data().arraydata()
labelW = labelWt[0].data().arraydata()
data = tex[0].np
idxW = idxWt[0].data().np
majorW = majorWt[0].data().np
labelW = labelWt[0].data().np

if size_min != 0 or depth_min != 0:
size = numpy.asarray([len(numpy.where(labelW == x)[0])
Expand All @@ -173,7 +170,7 @@ def watershedWithBasinsMerging(
# Blobs
blobindices = aimsalgo.blobsHeights(
mesh, tex[0].data(), labelWt[0].data())
blobheights = data[blobindices.arraydata()]
blobheights = data[blobindices.np]
blobdepths = data[idxW] - blobheights
print("criteria depth: ", blobdepths)

Expand Down
21 changes: 8 additions & 13 deletions pyaimsalgo/python/soma/wip/aimsalgo/foldsgraphthickness.py
Original file line number Diff line number Diff line change
Expand Up @@ -32,14 +32,11 @@
# The fact that you are presently reading this means that you have had
# knowledge of the CeCILL-B license and that you accept its terms.

from __future__ import print_function
from __future__ import absolute_import
from soma import aims
from soma import aimsalgo
import numpy


class FoldsGraphThickness(object):
class FoldsGraphThickness:
def __init__(self, fold_graph, lgw_vol, gm_wm_mesh, gm_lcr_mesh,
voronoi=None):
self.fold_graph = fold_graph
Expand Down Expand Up @@ -89,7 +86,7 @@ def printbucket(bck, vol, value):
f1.doit(seed, [-self.LCR_label, -self.GM_label], seed_label_list)
self.voronoi_vol = f1.voronoiVol()
print("Voronoi in White matter")
n = numpy.array(voronoi_vol, copy=False)
n = self.voronoi_vol.np
n[n == -1] = -100 # avoid value -1 which doesn't seem to work (!)
del n
f1.doit(self.voronoi_vol, [-100], seed_label_list)
Expand Down Expand Up @@ -127,23 +124,21 @@ def vorTexCreation(mesh):
def graphProcess(self):
voxel_size = self.lgw_vol.header()["voxel_size"]
voxel_volume = voxel_size[0]*voxel_size[1]*voxel_size[2]
coords = self.mid_interface.arraydata() != -1
data = self.mid_interface.arraydata()[coords]
coords = self.mid_interface.np != -1
data = self.mid_interface[coords]
self.fold_graph['thickness_mean'] = data.mean()
self.fold_graph['thickness_std'] = data.std()
coords = self.lgw_vol.arraydata() == self.GM_label
coords = self.lgw_vol.np == self.GM_label
GM_voxels = numpy.where(coords)[0].size
self.fold_graph['GM_volume'] = GM_voxels * voxel_volume
coords = self.lgw_vol.arraydata() == self.LCR_label
coords = self.lgw_vol.np == self.LCR_label
LCR_voxels = numpy.where(coords)[0].size
self.fold_graph['LCR_volume'] = LCR_voxels * voxel_volume
self.fold_graph['CSF_volume'] = LCR_voxels * voxel_volume

def vertexProcess(self):
print("vertexProcess")
voxel_size = self.lgw_vol.header()["voxel_size"]
voxel_volume = voxel_size[0]*voxel_size[1]*voxel_size[2]
coords = self.mid_interface.arraydata() != -1

def vertex_mesh(self, v, skel):
cut_mesh = aims.SurfaceManip.meshExtract(
Expand All @@ -168,8 +163,8 @@ def vertex_mesh(self, v, skel):
for v in self.fold_graph.vertices():
try:
skel = v['skeleton_label']
print("skel : " + str(skel) + " compteur : " + str(compteur) + '/'
+ str(n))
print("skel : " + str(skel) + " compteur : " + str(compteur)
+ '/' + str(n))
compteur = compteur + 1
vertex_mesh(self, v, skel)
except:
Expand Down
Loading

0 comments on commit da9f0dc

Please sign in to comment.