diff --git a/docs/preprocessing.rst b/docs/preprocessing.rst
index 6a97dc87..a7949810 100644
--- a/docs/preprocessing.rst
+++ b/docs/preprocessing.rst
@@ -384,12 +384,10 @@ Processing the *Subject Anatomical Reference* T1w or T2w images
                               template='MNI152NLin2009cAsym',
                               output_spaces=['T1w'],
                               output_resolution=1.25,
-                              skull_strip_template='OASIS',
                               force_spatial_normalization=True,
                               freesurfer=True,
                               longitudinal=False,
                               debug=False,
-                              hires=True,
                               num_t1w=1)
 
 As of version 0.18 QSIPrep has been changed to be very flexible with anatomical
diff --git a/qsiprep/cli/run.py b/qsiprep/cli/run.py
index 902f6675..5ca36878 100644
--- a/qsiprep/cli/run.py
+++ b/qsiprep/cli/run.py
@@ -328,6 +328,12 @@ def get_parser():
         choices=['MNI152NLin2009cAsym'],
         default='MNI152NLin2009cAsym',
         help='volume template space (default: MNI152NLin2009cAsym)')
+    g_conf.add_argument(
+        "--native-template",
+        required=False,
+        action='store',
+        help="A preprocessed native space image from fmriprep or previous qsiprep"
+        "that will be used instead of AC-PC to define subject native space.")
     g_conf.add_argument(
         '--output-resolution', '--output_resolution',
         action='store',
@@ -411,17 +417,6 @@ def get_parser():
     # ANTs options
     g_ants = parser.add_argument_group(
         'Specific options for ANTs registrations')
-    g_ants.add_argument(
-        '--skull-strip-template', '--skull_strip_template',
-        action='store',
-        default='OASIS',
-        choices=['OASIS', 'NKI'],
-        help='select ANTs skull-stripping template (default: OASIS)')
-    g_ants.add_argument(
-        '--skull-strip-fixed-seed', '--skull_strip_fixed_seed',
-        action='store_true',
-        help='do not use a random seed for skull-stripping - will ensure '
-        'run-to-run replicability when used with --omp-nthreads 1')
     g_ants.add_argument(
         '--skip-anat-based-spatial-normalization', '--skip_anat_based_spatial_normalization',
         action='store_true',
@@ -953,7 +948,6 @@ def build_qsiprep_workflow(opts, retval):
         work_dir=work_dir,
         output_dir=str(output_dir),
         ignore=opts.ignore,
-        hires=False,
         anatomical_contrast=opts.anat_modality,
         freesurfer=opts.do_reconall,
         bids_filters=opts.bids_filters,
@@ -975,11 +969,10 @@ def build_qsiprep_workflow(opts, retval):
         denoise_before_combining=not opts.denoise_after_combining,
         write_local_bvecs=opts.write_local_bvecs,
         omp_nthreads=omp_nthreads,
-        skull_strip_template=opts.skull_strip_template,
-        skull_strip_fixed_seed=opts.skull_strip_fixed_seed,
         force_spatial_normalization=force_spatial_normalization,
         output_resolution=opts.output_resolution,
         template=opts.anatomical_template,
+        native_template=opts.native_template,
         bids_dir=bids_dir,
         motion_corr_to=opts.b0_motion_corr_to,
         hmc_transform=opts.hmc_transform,
diff --git a/qsiprep/interfaces/anatomical.py b/qsiprep/interfaces/anatomical.py
index fdf4a523..376594ad 100644
--- a/qsiprep/interfaces/anatomical.py
+++ b/qsiprep/interfaces/anatomical.py
@@ -10,6 +10,7 @@
 """
 
 from pkg_resources import resource_filename as pkgr
+from pathlib import Path
 import os.path as op
 import numpy as np
 from nipype import logging
@@ -18,7 +19,7 @@
 from scipy.spatial import distance
 from scipy import ndimage
 from nipype.interfaces.base import (traits, TraitedSpec, BaseInterfaceInputSpec,
-                                    SimpleInterface, File)
+                                    SimpleInterface, File, isdefined)
 from nipype.utils.filemanip import fname_presuffix
 from dipy.segment.threshold import otsu
 import nilearn.image as nim
@@ -422,8 +423,7 @@ class _GetTemplateInputSpec(BaseInterfaceInputSpec):
         'MNI152NLin2009cAsym',
         usedefault=True,
         mandatory=True)
-    t1_file = File(exists=True)
-    t2_file = File(exists=True)
+    native_template_file = traits.Str()
     mask_file = File(exists=True)
     infant_mode = traits.Bool(False, usedefault=True)
     anatomical_contrast = traits.Enum("T1w", "T2w", "none")
@@ -434,6 +434,8 @@ class _GetTemplateOutputSpec(BaseInterfaceInputSpec):
     template_file = File(exists=True)
     template_brain_file = File(exists=True)
     template_mask_file = File(exists=True)
+    native_template_file = File(exists=True)
+    native_template_mask_file = File(exists=True)
 
 
 class GetTemplate(SimpleInterface):
@@ -471,4 +473,19 @@ def _run_interface(self, runtime):
         else:
             raise NotImplementedError("Arbitrary templates not available yet")
 
+        # By default set the native template file and mask to the template file and template mask
+        self._results["native_template_file"] = ref_img
+        self._results["native_template_mask_file"] = ref_img_mask
+        if isdefined(self.inputs.native_template_file):
+            native_template = Path(self.inputs.native_template_file)
+            if native_template.exists():
+                self._results["native_template_file"] = str(
+                    native_template.absolute())
+                self._results["native_template_brain_file"] = _get_corresponding_mask(
+                    self._results["native_template_file"])
+
         return runtime
+
+
+def _get_corresponding_mask(anat_file):
+    return anat_file
\ No newline at end of file
diff --git a/qsiprep/workflows/anatomical/volume.py b/qsiprep/workflows/anatomical/volume.py
index 952a1262..518e211c 100644
--- a/qsiprep/workflows/anatomical/volume.py
+++ b/qsiprep/workflows/anatomical/volume.py
@@ -49,7 +49,7 @@ class DerivativesDataSink(FDerivativesDataSink):
 FS_VERSION="7.3.1"
 
 #  pylint: disable=R0914
-def init_anat_preproc_wf(template, debug, dwi_only,
+def init_anat_preproc_wf(template, native_template, debug, dwi_only,
                          infant_mode, longitudinal, omp_nthreads,
                          output_dir, num_anat_images, output_resolution,
                          nonlinear_register_to_template,
@@ -76,12 +76,14 @@ def init_anat_preproc_wf(template, debug, dwi_only,
                                   output_dir='.',
                                   anatomical_contrast="T1w",
                                   template='MNI152NLin2009cAsym',
+                                  native_template=None,
                                   output_resolution=1.25,
                                   dwi_only=False,
                                   infant_mode=False,
                                   nonlinear_register_to_template=True,
                                   longitudinal=False,
                                   debug=False,
+                                  native_reference=None,
                                   num_anat_images=1)
 
     **Parameters**
@@ -98,7 +100,8 @@ def init_anat_preproc_wf(template, debug, dwi_only,
             sure to choose a robust option for ``interpolation`` to avoid ringing
             artifacts. One option is 'BSpline', which matches mrtrix.
         template : str
-            Name of template targeted by ``template`` output space
+            Name of template targeted by ``template`` output space.
+        native_template : str or None
         debug : bool
             Enable debugging outputs
         longitudinal : bool
@@ -110,6 +113,9 @@ def init_anat_preproc_wf(template, debug, dwi_only,
             Directory in which to save reportlets
         output_dir : str
             Directory in which to save derivatives
+        native_reference : str or None
+            Image containing the reference for AC-PC alignment. By default is
+            just the template image.
         name : str, optional
             Workflow name (default: anat_preproc_wf)
 
@@ -167,7 +173,8 @@ def init_anat_preproc_wf(template, debug, dwi_only,
     get_template_image = pe.Node(
         GetTemplate(template_name=template,
                     infant_mode=infant_mode,
-                    anatomical_contrast=anatomical_contrast),
+                    anatomical_contrast=anatomical_contrast,
+                    native_template_file=str(native_template)),
         name="get_template_image")
 
     # Create the output reference grid_image
@@ -250,6 +257,9 @@ def init_anat_preproc_wf(template, debug, dwi_only,
         has_rois=has_rois,
         name='anat_normalization_wf')
 
+    if native_template is not None:
+        LOGGER.info("An special reference image was used to define subject native space")
+
     # Resampling
     rigid_acpc_resample_brain = pe.Node(
         ants.ApplyTransforms(input_image_type=0,
@@ -344,6 +354,8 @@ def init_anat_preproc_wf(template, debug, dwi_only,
         (anat_reference_wf, anat_normalization_wf, [
             ('outputnode.bias_corrected', 'inputnode.anatomical_reference')]),
         (get_template_image, anat_normalization_wf, [
+            ('native_template_file', 'inputnode.native_template_image'),
+            ('native_template_mask_file', 'inputnode.native_template_mask'),
             ('template_file', 'inputnode.template_image'),
             ('template_mask_file', 'inputnode.template_mask')]),
         (anat_normalization_wf, outputnode, [
@@ -677,8 +689,6 @@ def init_anat_normalization_wf(sloppy, template_name, omp_nthreads,
     Parameters
         template_image: str
             Path to an image that will be used for Rigid ACPC align
-        skull_strip_template : str
-            Name of ANTs skull-stripping template ('OASIS' or 'NKI')
         debug : bool
             Enable debugging outputs
         omp_nthreads : int
@@ -705,6 +715,7 @@ def init_anat_normalization_wf(sloppy, template_name, omp_nthreads,
     inputnode = pe.Node(
         niu.IdentityInterface(
             fields=['template_image', 'template_mask',
+                    'native_template_image', 'native_template_mask',
                     'anatomical_reference', 'brain_mask', 'roi']),
         name='inputnode')
     outputnode = pe.Node(
@@ -739,8 +750,8 @@ def init_anat_normalization_wf(sloppy, template_name, omp_nthreads,
 
     workflow.connect([
         (inputnode, acpc_reg, [
-            ('template_image', 'reference_image'),
-            ('template_mask', 'reference_mask'),
+            ('native_template_image', 'reference_image'),
+            ('native_template_mask', 'reference_mask'),
             ('anatomical_reference', 'moving_image'),
             ('roi', 'lesion_mask'),
             ('brain_mask', 'moving_mask')]),
@@ -795,10 +806,10 @@ def init_anat_normalization_wf(sloppy, template_name, omp_nthreads,
             ('template_image', 'reference_image'),
             ('template_mask', 'reference_mask')]),
         (inputnode, rigid_acpc_resample_mask, [
-            ('template_image', 'reference_image'),
+            ('native_template_image', 'reference_image'),
             ('brain_mask', 'input_image')]),
         (inputnode, rigid_acpc_resample_anat, [
-            ('template_image', 'reference_image'),
+            ('native_template_image', 'reference_image'),
             ('anatomical_reference', 'input_image')]),
         (extract_rigid_transform, rigid_acpc_resample_anat, [
             ('rigid_transform', 'transforms')]),
diff --git a/qsiprep/workflows/base.py b/qsiprep/workflows/base.py
index 1ca3948e..c55ce2b4 100644
--- a/qsiprep/workflows/base.py
+++ b/qsiprep/workflows/base.py
@@ -43,15 +43,15 @@
 
 def init_qsiprep_wf(
         subject_list, run_uuid, work_dir, output_dir, bids_dir, ignore, debug,
-        low_mem, anat_only, dwi_only, longitudinal, b0_threshold, hires,
+        low_mem, anat_only, dwi_only, longitudinal, b0_threshold,
         anatomical_contrast, denoise_before_combining, dwi_denoise_window,
         denoise_method, unringing_method, b1_biascorrect_stage,
         no_b0_harmonization, output_resolution, infant_mode, combine_all_dwis,
         distortion_group_merge, pepolar_method, omp_nthreads, bids_filters,
-        force_spatial_normalization, skull_strip_template,
-        skull_strip_fixed_seed, freesurfer, hmc_model, impute_slice_threshold,
+        force_spatial_normalization,
+        freesurfer, hmc_model, impute_slice_threshold,
         hmc_transform, shoreline_iters, eddy_config, write_local_bvecs,
-        template, motion_corr_to, b0_to_t1w_transform,
+        template, native_template, motion_corr_to, b0_to_t1w_transform,
         intramodal_template_iters, intramodal_template_transform,
         prefer_dedicated_fmaps, fmap_bspline, fmap_demean, use_syn, force_syn,
         raw_image_sdc):
@@ -81,7 +81,6 @@ def init_qsiprep_wf(
                               longitudinal=False,
                               b0_threshold=100,
                               freesurfer=False,
-                              hires=False,
                               denoise_before_combining=True,
                               dwi_denoise_window=7,
                               denoise_method='patch2self',
@@ -94,9 +93,8 @@ def init_qsiprep_wf(
                               omp_nthreads=1,
                               output_resolution=2.0,
                               hmc_model='3dSHORE',
-                              skull_strip_template='OASIS',
-                              skull_strip_fixed_seed=False,
                               template='MNI152NLin2009cAsym',
+                              native_template=None,
                               motion_corr_to='iterative',
                               b0_to_t1w_transform='Rigid',
                               intramodal_template_iters=0,
@@ -169,15 +167,8 @@ def init_qsiprep_wf(
             Either 'DRBUDDI' or 'TOPUP'. The method for SDC when EPI fieldmaps are used.
         omp_nthreads : int
             Maximum number of threads an individual process may use
-        skull_strip_template : str
-            Name of ANTs skull-stripping template ('OASIS' or 'NKI')
-        skull_strip_fixed_seed : bool
-            Do not use a random seed for skull-stripping - will ensure
-            run-to-run replicability when used with --omp-nthreads 1
         freesurfer : bool
             Enable FreeSurfer surface reconstruction (may increase runtime)
-        hires : bool
-            Enable sub-millimeter preprocessing in FreeSurfer
         template : str
             Name of template targeted by ``template`` output space
         motion_corr_to : str
@@ -250,14 +241,12 @@ def init_qsiprep_wf(
             infant_mode=infant_mode,
             b0_threshold=b0_threshold,
             freesurfer=freesurfer,
-            hires=hires,
             combine_all_dwis=combine_all_dwis,
             distortion_group_merge=distortion_group_merge,
             pepolar_method=pepolar_method,
             omp_nthreads=omp_nthreads,
-            skull_strip_template=skull_strip_template,
-            skull_strip_fixed_seed=skull_strip_fixed_seed,
             template=template,
+            native_template=native_template,
             prefer_dedicated_fmaps=prefer_dedicated_fmaps,
             motion_corr_to=motion_corr_to,
             b0_to_t1w_transform=b0_to_t1w_transform,
@@ -290,9 +279,9 @@ def init_single_subject_wf(
         b0_threshold, denoise_before_combining, bids_filters, anatomical_contrast,
         dwi_denoise_window, denoise_method, unringing_method, b1_biascorrect_stage,
         no_b0_harmonization, infant_mode, combine_all_dwis, raw_image_sdc,
-        distortion_group_merge, pepolar_method, omp_nthreads, skull_strip_template,
-        force_spatial_normalization, skull_strip_fixed_seed, freesurfer, hires,
-        template, output_resolution, prefer_dedicated_fmaps,
+        distortion_group_merge, pepolar_method, omp_nthreads,
+        force_spatial_normalization, freesurfer,
+        template, native_template, output_resolution, prefer_dedicated_fmaps,
         motion_corr_to, b0_to_t1w_transform, intramodal_template_iters,
         intramodal_template_transform, hmc_model, hmc_transform,
         shoreline_iters, eddy_config, impute_slice_threshold, fmap_bspline,
@@ -336,15 +325,13 @@ def init_single_subject_wf(
             longitudinal=False,
             b0_threshold=100,
             freesurfer=False,
-            hires=False,
             force_spatial_normalization=True,
             combine_all_dwis=True,
             distortion_group_merge='none',
             pepolar_method='TOPUP',
             omp_nthreads=1,
-            skull_strip_template='OASIS',
-            skull_strip_fixed_seed=False,
             template='MNI152NLin2009cAsym',
+            native_template=None,
             prefer_dedicated_fmaps=False,
             motion_corr_to='iterative',
             b0_to_t1w_transform='Rigid',
@@ -407,15 +394,8 @@ def init_single_subject_wf(
             Either 'DRBUDDI' or 'TOPUP'. The method for SDC when EPI fieldmaps are used.
         omp_nthreads : int
             Maximum number of threads an individual process may use
-        skull_strip_template : str
-            Name of ANTs skull-stripping template ('OASIS' or 'NKI')
-        skull_strip_fixed_seed : bool
-            Do not use a random seed for skull-stripping - will ensure
-            run-to-run replicability when used with --omp-nthreads 1
         freesurfer : bool
             Enable FreeSurfer surface reconstruction (may increase runtime)
-        hires : bool
-            Enable sub-millimeter preprocessing in FreeSurfer
         reportlets_dir : str
             Directory in which to save reportlets
         output_dir : str
@@ -424,6 +404,8 @@ def init_single_subject_wf(
             Root directory of BIDS dataset
         template : str
             Name of template targeted by ``template`` output space
+        native_template : str or None
+            Either a path to an existing image or None
         hmc_model : 'none', '3dSHORE' or 'eddy'
             Model used to generate target images for head motion correction. If 'none'
             the transform from the nearest b0 will be used.
@@ -575,6 +557,7 @@ def init_single_subject_wf(
     info_modality = "dwi" if dwi_only else anatomical_contrast.lower()
     anat_preproc_wf = init_anat_preproc_wf(
         template=template,
+        native_template=native_template,
         debug=debug,
         dwi_only=dwi_only,
         infant_mode=infant_mode,
diff --git a/tests/get_data.py b/tests/get_data.py
index 4fcd114d..003db6d2 100644
--- a/tests/get_data.py
+++ b/tests/get_data.py
@@ -51,7 +51,6 @@ def get_default_cli_args():
         work_dir='',
         ignore=[],
         recon_only=False,
-        hires=False,
         freesurfer=False,
         do_reconall=False,
         debug=True,
@@ -63,8 +62,6 @@ def get_default_cli_args():
         denoise_before_combining=True,
         write_local_bvecs=False,
         omp_nthreads=1,
-        skull_strip_template="OASIS",
-        skull_strip_fixed_seed=False,
         force_spatial_normalization=False,
         output_resolution=5,
         template="MNI152NLin2009cAsym",
diff --git a/tests/smoke_tests.py b/tests/smoke_tests.py
index ec63672d..737c514f 100644
--- a/tests/smoke_tests.py
+++ b/tests/smoke_tests.py
@@ -22,14 +22,11 @@ def test_preproc_pepolar_sdc(tmp_path):
                          anat_only=False,
                          longitudinal=False,
                          freesurfer=False,
-                         hires=False,
                          force_spatial_normalization=False,
                          denoise_before_combining=True,
                          dwi_denoise_window=7,
                          combine_all_dwis=True,
                          omp_nthreads=1,
-                         skull_strip_template='OASIS',
-                         skull_strip_fixed_seed=False,
                          output_spaces=['T1w', 'template'],
                          template='MNI152NLin2009cAsym',
                          output_resolution=1.25,
@@ -66,14 +63,11 @@ def test_preproc_syn_sdc(tmp_path):
                          anat_only=False,
                          longitudinal=False,
                          freesurfer=False,
-                         hires=False,
                          denoise_before_combining=True,
                          dwi_denoise_window=7,
                          combine_all_dwis=True,
                          omp_nthreads=1,
-                         skull_strip_template='OASIS',
                          force_spatial_normalization=False,
-                         skull_strip_fixed_seed=False,
                          output_spaces=['T1w', 'template'],
                          template='MNI152NLin2009cAsym',
                          output_resolution=1.25,
@@ -109,13 +103,10 @@ def test_preproc_nonehmc_sdc(tmp_path):
                          anat_only=False,
                          longitudinal=False,
                          freesurfer=False,
-                         hires=False,
                          denoise_before_combining=True,
                          dwi_denoise_window=7,
                          combine_all_dwis=True,
                          omp_nthreads=1,
-                         skull_strip_template='OASIS',
-                         skull_strip_fixed_seed=False,
                          output_spaces=['T1w', 'template'],
                          force_spatial_normalization=False,
                          template='MNI152NLin2009cAsym',
@@ -153,14 +144,11 @@ def test_preproc_buds(tmp_path):
                          anat_only=False,
                          longitudinal=False,
                          freesurfer=False,
-                         hires=False,
                          denoise_before_combining=True,
                          force_spatial_normalization=False,
                          dwi_denoise_window=7,
                          combine_all_dwis=True,
                          omp_nthreads=1,
-                         skull_strip_template='OASIS',
-                         skull_strip_fixed_seed=False,
                          output_spaces=['T1w', 'template'],
                          template='MNI152NLin2009cAsym',
                          output_resolution=1.25,