With Rick's help, we were able to get around this snag.

What was going on was that I had been inserting blank slices as place holders which ended up moving my origin - this was the source of my coregistration issues.

I was doing this because 3dNLfim (which I was implementing slice by slice) does not output anything if the mask being used is all zeros (which some of my slices were). To get around this, I created a simple ROI of a line in the corner of the dataset that went through each slice and using 3dcalc, added this ROI to the mask I was already using.

While this means modeling dummy voxels, it's a neglible amount of extra processing time. Despite the fact that I am sure that more elegant solutions exists, this one was quick and dirty and seems to get the trick done.

Rick had suggested adding voxels to the empty slices using the following 3dcalc expression.

"So the expression was:

'isequal(i,0)*isequal(k,0)'

That would set one line to 1 in the case of slices
over the y direction. You might just add that to
an already existing mask:

3dcalc -a mask+orig -prefix mask_w_line \
-expr 'step(a + isequal(i,0)*isequal(k,0))'

So 'a' would be the old mask value, but then you would
also set any voxel along the line."