Python IOUtils.write_surface方法代码示例

# 需要导入模块: from tvb.recon.io.factory import IOUtils [as 别名]
# 或者: from tvb.recon.io.factory.IOUtils import write_surface [as 别名]
def test_write_write_brain_visa_surf():
    surface_path = get_data_file(
        TEST_FS_SUBJECT, TEST_SURFACE_FOLDER, "lh.pial")
    out_path = get_temporary_files_path("lh.pial.tri")

    surface = IOUtils.read_surface(surface_path, False)
    IOUtils.write_surface(out_path, surface)

    assert os.path.exists(out_path)

# 需要导入模块: from tvb.recon.io.factory import IOUtils [as 别名]
# 或者: from tvb.recon.io.factory.IOUtils import write_surface [as 别名]
def test_write_fs_surface():
    file_path = get_data_file(TEST_FS_SUBJECT, TEST_SURFACE_FOLDER, "lh.pial")
    original_surface = IOUtils.read_surface(file_path, False)
    triangles_number = len(original_surface.triangles)

    output_file_path = get_temporary_files_path("lh-test.pial")
    IOUtils.write_surface(output_file_path, original_surface)

    new_surface = IOUtils.read_surface(output_file_path, False)
    assert triangles_number == len(new_surface.triangles) == 327680

# 需要导入模块: from tvb.recon.io.factory import IOUtils [as 别名]
# 或者: from tvb.recon.io.factory.IOUtils import write_surface [as 别名]
    def aseg_surf_conc_annot(self, surf_path: str, out_surf_path: str, annot_path: str,
                             label_indices: Union[numpy.ndarray, list], lut_path: Optional[str]=None) -> Surface:
        """
        Concatenate surfaces of one specific label of interest each, to create a single annotated surface.
        """

        lut_path = lut_path or default_lut_path()

        label_names, color_table = self.annotation_service.lut_to_annot_names_ctab(lut_path=lut_path,
                                                                                   labels=label_indices)
        label_indices = numpy.array(label_indices.split()).astype('i')

        #                  verts tri area_mask cras
        surfaces = []
        out_annotation = Annotation([], [], [])
        label_number = -1

        for label_index in label_indices:
            # TODO: This is hardcoded: /aseg-%06d here and also in pegasus dax generator
            this_surf_path = surf_path + "/aseg-%06d" % int(label_index)

            if os.path.exists(this_surf_path):
                ind_l, = numpy.where(label_indices == label_index)
                out_annotation.add_region_names_and_colors(
                    label_names[int(ind_l)],
                    color_table[ind_l, :])
                label_number += 1
                surfaces.append(IOUtils.read_surface(this_surf_path, False))
                out_annotation.add_region_mapping(
                    label_number * numpy.ones((surfaces[-1].n_vertices,), dtype='int64'))
        out_surface = self.merge_surfaces(surfaces)
        # out_annotation.regions_color_table = numpy.squeeze(numpy.array(out_annotation.regions_color_table).astype('i'))

        IOUtils.write_surface(out_surf_path, out_surface)
        IOUtils.write_annotation(annot_path, out_annotation)

        return out_surface

# 需要导入模块: from tvb.recon.io.factory import IOUtils [as 别名]
# 或者: from tvb.recon.io.factory.IOUtils import write_surface [as 别名]
def compute_region_details(atlas_suffix: AtlasSuffix, fs_color_lut: os.PathLike, t1: os.PathLike, lh_cort: os.PathLike,
                           rh_cort: os.PathLike, lh_cort_annot: os.PathLike, rh_cort_annot: os.PathLike,
                           lh_subcort: os.PathLike, rh_subcort: os.PathLike, lh_subcort_annot: os.PathLike,
                           rh_subcort_annot: os.PathLike):
    annot_cort_lh = IOUtils.read_annotation(lh_cort_annot)
    annot_cort_rh = IOUtils.read_annotation(rh_cort_annot)

    annot_subcort_lh = IOUtils.read_annotation(lh_subcort_annot)
    annot_subcort_rh = IOUtils.read_annotation(rh_subcort_annot)

    mapping = MappingService(atlas_suffix, annot_cort_lh, annot_cort_rh, annot_subcort_lh, annot_subcort_rh)
    mapping.generate_region_mapping_for_cort_annot(annot_cort_lh, annot_cort_rh)
    mapping.generate_region_mapping_for_subcort_annot(annot_subcort_lh, annot_subcort_rh)

    surface_service = SurfaceService()

    surf_cort_lh = IOUtils.read_surface(lh_cort, False)
    surf_cort_rh = IOUtils.read_surface(rh_cort, False)

    full_cort_surface = surface_service.merge_surfaces([surf_cort_lh, surf_cort_rh])

    surf_subcort_lh = IOUtils.read_surface(lh_subcort, False)
    surf_subcort_rh = IOUtils.read_surface(rh_subcort, False)

    full_subcort_surface = surface_service.merge_surfaces([surf_subcort_lh, surf_subcort_rh])

    genericIO.write_list_to_txt_file(mapping.cort_region_mapping, AsegFiles.RM_CORT_TXT.value.replace("%s", atlas_suffix))
    genericIO.write_list_to_txt_file(mapping.subcort_region_mapping,
                                     AsegFiles.RM_SUBCORT_TXT.value.replace("%s", atlas_suffix))

    vox2ras_file = "vox2ras.txt"
    subprocess.call(["mri_info", "--vox2ras", t1, "--o", vox2ras_file])

    surf_subcort_filename = "surface_subcort.zip"
    IOUtils.write_surface(surf_subcort_filename, full_subcort_surface)

    surf_cort_filename = "surface_cort.zip"
    IOUtils.write_surface(surf_cort_filename, full_cort_surface)

    os.remove(vox2ras_file)

    cort_subcort_full_surf = surface_service.merge_surfaces([full_cort_surface, full_subcort_surface])
    cort_subcort_full_region_mapping = mapping.cort_region_mapping + mapping.subcort_region_mapping

    dict_fs_custom = mapping.get_mapping_for_connectome_generation()
    genericIO.write_dict_to_txt_file(dict_fs_custom, AsegFiles.FS_CUSTOM_TXT.value.replace("%s", atlas_suffix))

    region_areas = surface_service.compute_areas_for_regions(mapping.get_all_regions(), cort_subcort_full_surf,
                                                             cort_subcort_full_region_mapping)
    genericIO.write_list_to_txt_file(region_areas, AsegFiles.AREAS_TXT.value.replace("%s", atlas_suffix))

    region_centers = surface_service.compute_centers_for_regions(mapping.get_all_regions(), cort_subcort_full_surf,
                                                                 cort_subcort_full_region_mapping)
    cort_subcort_lut = mapping.get_entire_lut()
    region_names = list(cort_subcort_lut.values())

    with open(AsegFiles.CENTERS_TXT.value.replace("%s", atlas_suffix), "w") as f:
        for idx, (val_x, val_y, val_z) in enumerate(region_centers):
            f.write("%s %.2f %.2f %.2f\n" % (region_names[idx], val_x, val_y, val_z))

    region_orientations = surface_service.compute_orientations_for_regions(mapping.get_all_regions(),
                                                                           cort_subcort_full_surf,
                                                                           cort_subcort_full_region_mapping)

    lh_region_centers = surface_service.compute_centers_for_regions(mapping.get_lh_regions(), surf_cort_lh,
                                                                    mapping.lh_region_mapping)
    lh_region_orientations = surface_service.compute_orientations_for_regions(mapping.get_lh_regions(), surf_cort_lh,
                                                                              mapping.lh_region_mapping)
    with open(AsegFiles.LH_DIPOLES_TXT.value.replace("%s", atlas_suffix), "w") as f:
        for idx, (val_x, val_y, val_z) in enumerate(lh_region_centers):
            f.write("%.2f %.2f %.2f %.2f %.2f %.2f\n" % (
                val_x, val_y, val_z, lh_region_orientations[idx][0], lh_region_orientations[idx][1],
                lh_region_orientations[idx][2]))

    rh_region_centers = surface_service.compute_centers_for_regions(mapping.get_rh_regions(), surf_cort_rh,
                                                                    mapping.rh_region_mapping)
    rh_region_orientations = surface_service.compute_orientations_for_regions(mapping.get_rh_regions(), surf_cort_rh,
                                                                              mapping.rh_region_mapping)
    with open(AsegFiles.RH_DIPOLES_TXT.value.replace("%s", atlas_suffix), "w") as f:
        for idx, (val_x, val_y, val_z) in enumerate(rh_region_centers):
            f.write("%.2f %.2f %.2f %.2f %.2f %.2f\n" % (
                val_x, val_y, val_z, rh_region_orientations[idx][0], rh_region_orientations[idx][1],
                rh_region_orientations[idx][2]))

    numpy.savetxt(AsegFiles.ORIENTATIONS_TXT.value.replace("%s", atlas_suffix), region_orientations, fmt='%.2f %.2f %.2f')

    annotation_service = AnnotationService()
    lut_dict, _, _ = annotation_service.read_lut(fs_color_lut, "name")
    rm_index_dict = mapping.get_mapping_for_aparc_aseg(lut_dict)
    genericIO.write_dict_to_txt_file(rm_index_dict, AsegFiles.RM_TO_APARC_ASEG_TXT.value.replace("%s", atlas_suffix))

    genericIO.write_list_to_txt_file(mapping.is_cortical_region_mapping(),
                                     AsegFiles.CORTICAL_TXT.value.replace("%s", atlas_suffix))

# 需要导入模块: from tvb.recon.io.factory import IOUtils [as 别名]
# 或者: from tvb.recon.io.factory.IOUtils import write_surface [as 别名]
def generate_surface_zip(in_file, out_file):
    surface = IOUtils.read_surface(in_file, False)
    IOUtils.write_surface(out_file, surface)

# 需要导入模块: from tvb.recon.io.factory import IOUtils [as 别名]
# 或者: from tvb.recon.io.factory.IOUtils import write_surface [as 别名]
def merge_surfs(surf_lh, surf_rh, out_surf_path):
    s_lh = IOUtils.read_surface(surf_lh, False)
    s_rh = IOUtils.read_surface(surf_rh, False)
    surf = surfaceService.merge_surfaces([s_lh, s_rh])
    IOUtils.write_surface(out_surf_path, surf)

# 需要导入模块: from tvb.recon.io.factory import IOUtils [as 别名]
# 或者: from tvb.recon.io.factory.IOUtils import write_surface [as 别名]
 def convert_fs_to_brain_visa(self, in_surf_path: str, out_surf_path: Optional[str]=None):
     surface = IOUtils.read_surface(in_surf_path, False)
     if out_surf_path is None:
         out_surf_path = in_surf_path + '.tri'
     IOUtils.write_surface(out_surf_path, surface)

# 需要导入模块: from tvb.recon.io.factory import IOUtils [as 别名]
# 或者: from tvb.recon.io.factory.IOUtils import write_surface [as 别名]

#.........这里部分代码省略.........
        # Initialize the output mask:
        verts_out_mask = numpy.repeat([False], surface.vertices.shape[0])
        for label_index in range(len(region_mapping_indexes)):

            self.logger.info("%s", add_string +
                             annotation.region_names[label_index])

            # Get the indexes of the vertices corresponding to this label:
            verts_indices_of_label, = numpy.where(
                annotation.region_mapping[:] == region_mapping_indexes[label_index])
            verts_indices_of_label_size = verts_indices_of_label.size
            if verts_indices_of_label_size == 0:
                continue

            # Add any additional labels
            all_labels = [labels[label_index]] + add_lbl

            # get the vertices for current label and add cras to take them to
            # scanner ras
            verts_of_label = surface.vertices[verts_indices_of_label, :]
            verts_of_label += numpy.repeat(numpy.expand_dims(
                cras, 1).T, verts_indices_of_label_size, axis=0)

            # Compute the nearest voxel coordinates using the affine transform
            ijk = numpy.round(
                ras2vox_affine_matrix.dot(numpy.c_[verts_of_label, numpy.ones(verts_indices_of_label_size)].T)[:3].T) \
                .astype('i')

            # Get the labels of these voxels:
            surf_vxls = volume.data[ijk[:, 0], ijk[:, 1], ijk[:, 2]]

            # Vertex mask to keep: those that correspond to voxels of one of
            # the target labels
            # surf_vxls==lbl if only one target label
            verts_keep, = numpy.where(numpy.in1d(surf_vxls, all_labels))
            verts_out_mask[verts_indices_of_label[verts_keep]] = True

            if vertex_neighbourhood > 0:
                # These are now the remaining indexes to be checked for
                # neighboring voxels
                verts_indices_of_label = numpy.delete(
                    verts_indices_of_label, verts_keep)
                ijk = numpy.delete(ijk, verts_keep, axis=0)

                for vertex_index in range(verts_indices_of_label.size):
                    # Generate the specific grid centered at the voxel ijk
                    ijk_grid = grid + \
                               numpy.tile(ijk[vertex_index, :], (n_grid, 1))

                    # Remove voxels outside the volume
                    indexes_within_limits = numpy.all([(ijk_grid[:, 0] >= 0), (ijk_grid[:, 0] < volume.dimensions[0]),
                                                       (ijk_grid[:, 1] >= 0), (ijk_grid[
                                                                               :, 1] < volume.dimensions[1]),
                                                       (ijk_grid[:, 2] >= 0), (ijk_grid[:, 2] < volume.dimensions[2])],
                                                      axis=0)
                    ijk_grid = ijk_grid[indexes_within_limits, :]
                    surf_vxls = volume.data[
                        ijk_grid[:, 0], ijk_grid[:, 1], ijk_grid[:, 2]]

                    # If any of the neighbors is of the target labels include
                    # the current vertex
                    # surf_vxls==lbl if only one target label
                    if numpy.any(numpy.in1d(surf_vxls, all_labels)):
                        verts_out_mask[
                            verts_indices_of_label[vertex_index]] = True

        # Vertex indexes and vertices to keep:
        verts_out_indices, = numpy.where(verts_out_mask)
        verts_out = surface.vertices[verts_out_indices]

        # TODO maybe: make sure that all voxels of this label correspond to at least one vertex.
        # Create a similar mask for faces by picking only triangles of which
        # all 3 vertices are included
        face_out_mask = numpy.c_[
            verts_out_mask[surface.triangles[:, 0]], verts_out_mask[surface.triangles[:, 1]], verts_out_mask[
                surface.triangles[:, 2]]].all(axis=1)
        faces_out = surface.triangles[face_out_mask]

        # The old vertices' indexes of faces have to be transformed to the new
        # vrtx_out_inds:
        for iF in range(faces_out.shape[0]):
            for vertex_index in range(3):
                faces_out[iF, vertex_index], = numpy.where(
                    faces_out[iF, vertex_index] == verts_out_indices)

        surface.vertices = verts_out
        surface.triangles = faces_out

        # Write the output surfaces and annotations to files. Also write files
        # with the indexes of vertices to keep.
        IOUtils.write_surface(out_surf_path, surface)

        annotation.set_region_mapping(
            annotation.get_region_mapping_by_indices([verts_out_indices]))
        IOUtils.write_annotation(out_surf_path + ".annot", annotation)

        numpy.save(out_surf_path + "-idx.npy", verts_out_indices)
        numpy.savetxt(out_surf_path + "-idx.txt", verts_out_indices, fmt='%d')

        return surface, annotation