An important point to remember is a half-voxel shift (the maximum interpolation) will introduce just as much interpolation as a seven and a half voxel shift, so you may not need to worry that much about differences between runs. The interpolation introduced by motion correction (volume registration) can very well be this large within a run.
If it is still a concern, you could use @align_centers to match the center of mass across datasets to some reference as a preprocessing step. This does not actually change any data at all; it only moves the origin in the header, so there is no interpolation at all. Then processing with @align_epi_anat.py using those datasets would work.
Other options could be to limit alignment to Nearest Neighbor interpolation, but that would effectively limit alignment to voxel increments. You could align anatomical data to each of these runs separately to find anatomical references.
If you will be transforming the data to match a template (Talairach) to get anatomical landmarks, the shifts will be larger and between run shifts will be less important. In that case, the align_epi_anat.py with a -tlrc_apar option can transform the epi and child_epi datasets to match the talaraich transformed anatomical dataset.