Hyperslabs part 2
editOnto more hyperslabs: we will once again do the same trick of adding one to every voxel. We'll do it slightly differently: this time we will get and set the hyperslab at every slice, using much less memory. And, for even more mincy goodness, we will use slice scaling for the output.
#include <minc2.h>
#include <stdio.h>
#include <stdlib.h>
#include <assert.h>
#include <string.h>
int main(int argc, char **argv) {
mihandle_t minc_volume, new_volume;
midimhandle_t dimensions[3], *dimensions_new;
double min, max;
int result, i, slice;
unsigned long start[3], count[3];
unsigned int sizes[3];
double *slab;
if (argc != 3) {
fprintf(stderr, "Usage: input.mnc output.mnc\n");
}
/* open the volume - first command line argument */
result = miopen_volume(argv[1], MI2_OPEN_READ, &minc_volume);
/* check for error on opening */
if (result != MI_NOERROR) {
fprintf(stderr, "Error opening input file: %d.\n", result);
return(1);
}
/* get the dimension sizes */
miget_volume_dimensions(minc_volume, MI_DIMCLASS_SPATIAL,
MI_DIMATTR_ALL, MI_DIMORDER_FILE,
3, dimensions);
miget_dimension_sizes(dimensions, 3, sizes);
printf("Volume sizes: %u %u %u\n", sizes[0], sizes[1], sizes[2]);
/* allocate new dimensions */
dimensions_new = malloc(sizeof(midimhandle_t) * 3);
/* copy the dimensions */
for(i=0; i < 3; i++) {
micopy_dimension(dimensions[i], &dimensions_new[i]);
}
This code is unchanged from the last tutorial.
/* create the new volume */
if (micreate_volume(argv[2], 3, dimensions_new, MI_TYPE_UBYTE,
MI_CLASS_REAL, NULL, &new_volume) != MI_NOERROR) {
fprintf(stderr, "Error creating new volume\n");
return(1);
}
/* indicate that we will be using slice scaling */
miset_slice_scaling_flag(new_volume, TRUE);
/* set valid range */
miset_volume_valid_range(new_volume, 255, 0);
/* create the data for the new volume */
if (micreate_volume_image(new_volume) != MI_NOERROR) {
fprintf(stderr, "Error creating volume data\n");
return(1);
}
This time we create the volume before we start messing with the actual data. We also set slice scaling using miset_slice_scaling_flag - note that this has to come before micreate_volume_image. The rest is the same as before.
/* now go into a slice loop - get each slice as a hyperslab. *
* First allocate enough memory to hold one slice of data */
slab = malloc(sizeof(double) * sizes[1] * sizes[2]);
/* the start and counts */
start[0] = start[1] = start[2] = 0;
/* always one slice */
count[0] = 1;
count[1] = (unsigned long) sizes[1];
count[2] = (unsigned long) sizes[2];
/* start the loop itself */
for (slice=0; slice < sizes[0]; slice++) {
/* set the slice start */
start[0] = (unsigned long) slice;
if (miget_real_value_hyperslab(minc_volume, MI_TYPE_DOUBLE,
start, count, slab) != MI_NOERROR) {
fprintf(stderr, "Error getting hyperslab for slice %d\n", slice);
return(1);
}
/* set min and max to the first voxel of the slice plus 1
* (since we add one to everything anyway) */
min = slab[0]+1;
max = slab[0]+1;
/* loop over all voxels in slice */
for (i=0; i < sizes[1] * sizes[2]; i++) {
/* increment voxel by 1 */
slab[i] += 1;
/* check if min or max should be changed */
if (slab[i]+1 < min) {
min = slab[i]+1;
}
else if (slab[i]+1 > max) {
max = slab[i]+1;
}
}
/* set the slice scaling */
miset_slice_range(new_volume, start, 3, max, min);
/* put the hyperslab into the new volume */
if (miset_real_value_hyperslab(new_volume, MI_TYPE_DOUBLE,
start, count, slab) != MI_NOERROR) {
fprintf(stderr, "Error setting hyperslab\n");
return(1);
}
}
The loop over the data has become more complicated. We now allocate data for a hyperslab big enough to hold one slice, and we reuse the same array for every slice. We then get the hyperslab for every slice, do the same manipulation of the voxels while keeping track of slice specific min and max, and set the slice range. Then we put the modified hyperslab into the new volume, and move on to the next slice.
/* closes the volume and makes sure all data is written to file */
miclose_volume(minc_volume);
miclose_volume(new_volume);
/* free memory */
free(dimensions_new);
free(slab);
return(0);
}
And we finish up in the same way as before.
code
edit#include <minc2.h>
#include <stdio.h>
#include <stdlib.h>
#include <assert.h>
#include <string.h>
int main(int argc, char **argv) {
mihandle_t minc_volume, new_volume;
midimhandle_t dimensions[3], *dimensions_new;
double min, max;
int result, i, slice;
unsigned long start[3], count[3];
unsigned int sizes[3];
double *slab;
if (argc != 3) {
fprintf(stderr, "Usage: input.mnc output.mnc\n");
}
/* open the volume - first command line argument */
result = miopen_volume(argv[1], MI2_OPEN_READ, &minc_volume);
/* check for error on opening */
if (result != MI_NOERROR) {
fprintf(stderr, "Error opening input file: %d.\n", result);
return(1);
}
/* get the dimension sizes */
miget_volume_dimensions(minc_volume, MI_DIMCLASS_SPATIAL,
MI_DIMATTR_ALL, MI_DIMORDER_FILE,
3, dimensions);
miget_dimension_sizes(dimensions, 3, sizes);
printf("Volume sizes: %u %u %u\n", sizes[0], sizes[1], sizes[2]);
/* allocate new dimensions */
dimensions_new = malloc(sizeof(midimhandle_t) * 3);
/* copy the dimensions */
for(i=0; i < 3; i++) {
micopy_dimension(dimensions[i], &dimensions_new[i]);
}
/* create the new volume */
if (micreate_volume(argv[2], 3, dimensions_new, MI_TYPE_UBYTE,
MI_CLASS_REAL, NULL, &new_volume) != MI_NOERROR) {
fprintf(stderr, "Error creating new volume\n");
return(1);
}
/* indicate that we will be using slice scaling */
miset_slice_scaling_flag(new_volume, TRUE);
/* set valid range */
miset_volume_valid_range(new_volume, 255, 0);
/* create the data for the new volume */
if (micreate_volume_image(new_volume) != MI_NOERROR) {
fprintf(stderr, "Error creating volume data\n");
return(1);
}
/* now go into a slice loop - get each slice as a hyperslab. *
* First allocate enough memory to hold one slice of data */
slab = malloc(sizeof(double) * sizes[1] * sizes[2]);
/* the start and counts */
start[0] = start[1] = start[2] = 0;
/* always one slice */
count[0] = 1;
count[1] = (unsigned long) sizes[1];
count[2] = (unsigned long) sizes[2];
/* start the loop itself */
for (slice=0; slice < sizes[0]; slice++) {
/* set the slice start */
start[0] = (unsigned long) slice;
if (miget_real_value_hyperslab(minc_volume, MI_TYPE_DOUBLE,
start, count, slab) != MI_NOERROR) {
fprintf(stderr, "Error getting hyperslab for slice %d\n", slice);
return(1);
}
/* set min and max to the first voxel of the slice plus 1
* (since we add one to everything anyway) */
min = slab[0]+1;
max = slab[0]+1;
/* loop over all voxels in slice */
for (i=0; i < sizes[1] * sizes[2]; i++) {
/* increment voxel by 1 */
slab[i] += 1;
/* check if min or max should be changed */
if (slab[i]+1 < min) {
min = slab[i]+1;
}
else if (slab[i]+1 > max) {
max = slab[i]+1;
}
}
/* set the slice scaling */
miset_slice_range(new_volume, start, 3, max, min);
/* put the hyperslab into the new volume */
if (miset_real_value_hyperslab(new_volume, MI_TYPE_DOUBLE,
start, count, slab) != MI_NOERROR) {
fprintf(stderr, "Error setting hyperslab\n");
return(1);
}
}
/* closes the volume and makes sure all data is written to file */
miclose_volume(minc_volume);
miclose_volume(new_volume);
/* free memory */
free(dimensions_new);
free(slab);
return(0);
}