MINC/Tutorials/Programming04

Writing a minc volume to file, with an aside on slice scaling

The code in this tutorial reads a MINC file, adds 1 to every voxel, and then writes the modified volume back to file.

Two notes before getting to the code: (1) this will overwrite the file passed as the first argument, so be careful if you actually try the code, and (2) this is a slow way of accomplishing the task - see the next tutorial for a much more sensible way of accomplishing the same goal.

On to the code:

#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;
  midimhandle_t dimensions[3];
  double        voxel, max;
  int           result;
  unsigned long v[3];
  unsigned int  sizes[3];
  miboolean_t   slice_scaling;

  /* open the volume - first command line argument */
  result = miopen_volume(argv[1], MI2_OPEN_RDWR, &minc_volume);
  /* check for error on opening */
  if (result != MI_NOERROR) {
    fprintf(stderr, "Error opening input file: %d.\n", result);
    return(1);
  }

Nothing new up to this point - the variables are declared, the volume is read (though note the second argument to miopen_volume - it now allows for both reading and writing of this volume).

  /* check whether this volume uses slice scaling */
  miget_slice_scaling_flag(minc_volume, &slice_scaling);

Here the first new minc function crops up: miget_slice_scaling_flag. It performs a very simple operation, asking whether the volume uses slice scaling or global scaling, placing the output into the boolean type miboolean_t.

Scaling is one of the key concepts in MINC, and is introduced here: MINC/Tutorials/Tutorials01 (scroll down to the data-types section). Our program needs to know whether slice scaling is being used for the following reason: the range represented in the volume is usually set to follow the actual range of values present in the volume. This means that if we increment all voxels by 1, we have to make sure that the range of values allowed in the volume is updated as well. As you will see below, updating the range is different depending on whether the scaling is set for the global volume or independently for every slice.

On to more source code:

  /* get the dimension sizes */
  miget_volume_dimensions(minc_volume, MI_DIMCLASS_SPATIAL,
			  MI_DIMATTR_ALL, MI_DIMORDER_FILE,
			  3, dimensions);
  result = miget_dimension_sizes(dimensions, 3, sizes);
  printf("Volume sizes: %u %u %u\n", sizes[0], sizes[1], sizes[2]);

We obtain the dimension sizes once again ...

  /* global scaling: increment the maximum by one */
  if (slice_scaling == FALSE) {
    miget_volume_max(minc_volume, &max);
    miset_volume_max(minc_volume, max+1);
  }

... and then proceed to the first part of the scaling issue. If there is no slice scaling, we just need to change the global scaling, which is done by getting the volume maximum using miget_volume_max, and setting it to the previous max+1 using miset_volume_max. There are similar functions for volume_min, which we could reasonably have updated here as well.

  /* loop over all voxels in the volume, add 1 to each voxel */
  for (v[0]=0; v[0] < sizes[0]; v[0]++) {
    /* slice scaling: increment the maximum of each slice by one */
    if(slice_scaling == TRUE) {
      miget_slice_max(minc_volume, v, 3, &max);
      miset_slice_max(minc_volume, v, 3, max+1);
    }
    for (v[1]=0; v[1] < sizes[1]; v[1]++) {
      for (v[2]=0; v[2] < sizes[2]; v[2]++) {
	miget_real_value(minc_volume, v, 3, &voxel);
	voxel += 1;
	miset_real_value(minc_volume, v, 3, voxel);
      }
    }
  }

This is the core of the code - we loop over all voxels using dimension size information to know when to terminate each loop, get the real value at each voxel, increment it by one, and write it back to the volume using miset_real_value. The first dimension is the slice dimension - so if this volume uses slice scaling we get the slice maximum using miget_slice_max, increment it by one and write it back using miset_slice_max.

  /* closes the volume and makes sure all data is written to file */
  miclose_volume(minc_volume);

  return(0);
}

Once we've looped over all voxels we close the volume using miclose_volume. And that's it, we're done. The file read in as the first argument should now have all of its values incremented by one. Keep in mind that this was the slow way of doing it - see the next tutorial for a much more sensible version.

code

#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;
  midimhandle_t dimensions[3];
  double        voxel, max;
  int           result;
  unsigned long v[3];
  unsigned int  sizes[3];
  miboolean_t   slice_scaling;

  /* open the volume - first command line argument */
  result = miopen_volume(argv[1], MI2_OPEN_RDWR, &minc_volume);
  /* check for error on opening */
  if (result != MI_NOERROR) {
    fprintf(stderr, "Error opening input file: %d.\n", result);
    return(1);
  }

  /* check whether this volume uses slice scaling */
  miget_slice_scaling_flag(minc_volume, &slice_scaling);

  /* get the dimension sizes */
  miget_volume_dimensions(minc_volume, MI_DIMCLASS_SPATIAL,
			  MI_DIMATTR_ALL, MI_DIMORDER_FILE,
			  3, dimensions);
  result = miget_dimension_sizes(dimensions, 3, sizes);
  printf("Volume sizes: %u %u %u\n", sizes[0], sizes[1], sizes[2]);

  /* global scaling: increment the maximum by one */
  if (slice_scaling == FALSE) {
    miget_volume_max(minc_volume, &max);
    miset_volume_max(minc_volume, max+1);
  }
  /* loop over all voxels in the volume, add 1 to each voxel */
  for (v[0]=0; v[0] < sizes[0]; v[0]++) {
    /* slice scaling: increment the maximum of each slice by one */
    if(slice_scaling == TRUE) {
      miget_slice_max(minc_volume, v, 3, &max);
      miset_slice_max(minc_volume, v, 3, max+1);
    }
    for (v[1]=0; v[1] < sizes[1]; v[1]++) {
      for (v[2]=0; v[2] < sizes[2]; v[2]++) {
	miget_real_value(minc_volume, v, 3, &voxel);
	voxel += 1;
	miset_real_value(minc_volume, v, 3, voxel);
      }
    }
  }
  
  /* closes the volume and makes sure all data is written to file */
  miclose_volume(minc_volume);

  return(0);
}