def saveinr(vol, h, fname):
    """
    Save a volume (described as a numpy array) to INR format.
    Code inspired by iso2mesh (http://iso2mesh.sf.net) by Q. Fang
    INPUTS:
    - vol: volume as numpy array
    - h: voxel sizes as list or numpy array
    - fname: filename for saving the inr file
    """
    fid = open(fname, 'wb')
    vol_type = vol.dtype
    if vol_type == 'uint8':
        btype = 'unsigned fixed'
        bitlen = 8
    elif vol_type == 'uint16':
        btype = 'unsigned fixed'
        bitlen = 16
    elif vol_type == 'float32':
        btype = 'float'
        bitlen = 32
    elif vol_type == 'float64':
        btype = 'float'
        bitlen = 64
    else:
        raise vol_type
    header = (
        '#INRIMAGE-4#{8:s}\nXDIM={0:d}\nYDIM={1:d}\nZDIM={2:d}\nVDIM=1\nTYPE={3:s}\n'
         + """PIXSIZE={4:d} bits
CPU=decm
VX={5:f}
VY={6:f}
VZ={7:f}
"""
        ).format(*vol.shape, btype, bitlen, h[0], h[1], h[2], '{')
    header = header + '\n' * (256 - 4 - len(header)) + '##}\n'
    fid.write(header.encode('ascii'))
    fid.write(vol.tobytes(order='F'))
