How to use the niworkflows.interfaces.ants.ImageMath function in niworkflows

'std2anat_xfm',
        'std_dseg',
        'std_mask',
        'std_tpms',
        'template',
        'template_spec',
    ]
    poutputnode = pe.Node(niu.IdentityInterface(fields=out_fields), name='poutputnode')

    split_desc = pe.Node(TemplateDesc(), run_without_submitting=True, name='split_desc')

    tf_select = pe.Node(TemplateFlowSelect(resolution=1 + debug),
                        name='tf_select', run_without_submitting=True)

    # With the improvements from poldracklab/niworkflows#342 this truncation is now necessary
    trunc_mov = pe.Node(ImageMath(operation='TruncateImageIntensity', op2='0.01 0.999 256'),
                        name='trunc_mov')

    registration = pe.Node(RobustMNINormalization(
        float=True, flavor=['precise', 'testing'][debug],
    ), name='registration', n_procs=omp_nthreads, mem_gb=2)

    # Resample T1w-space inputs
    tpl_moving = pe.Node(ApplyTransforms(
        dimension=3, default_value=0, float=True,
        interpolation='LanczosWindowedSinc'), name='tpl_moving')

    std_mask = pe.Node(ApplyTransforms(interpolation='MultiLabel'), name='std_mask')
    std_dseg = pe.Node(ApplyTransforms(interpolation='MultiLabel'), name='std_dseg')

    std_tpms = pe.MapNode(ApplyTransforms(dimension=3, default_value=0, float=True,
                                          interpolation='Gaussian'),

# De-pad
    depad_mask = pe.Node(
        ImageMath(operation="PadImage", op2="-%d" % padding), name="23_depad_mask"
    )
    depad_segm = pe.Node(
        ImageMath(operation="PadImage", op2="-%d" % padding), name="24_depad_segm"
    )
    depad_gm = pe.Node(
        ImageMath(operation="PadImage", op2="-%d" % padding), name="25_depad_gm"
    )
    depad_wm = pe.Node(
        ImageMath(operation="PadImage", op2="-%d" % padding), name="26_depad_wm"
    )
    depad_csf = pe.Node(
        ImageMath(operation="PadImage", op2="-%d" % padding), name="27_depad_csf"
    )

    msk_conform = pe.Node(niu.Function(function=_conform_mask), name="msk_conform")
    merge_tpms = pe.Node(niu.Merge(in_segmentation_model[0]), name="merge_tpms")
    # fmt: off
    wf.connect([
        (inputnode, copy_xform, [(("in_files", _pop), "hdr_file")]),
        (inputnode, pad_mask, [("in_mask", "op1")]),
        (inputnode, atropos, [
            ("in_files", "intensity_images"),
            ("in_mask_dilated", "mask_image"),
        ]),
        (inputnode, msk_conform, [(("in_files", _pop), "in_reference")]),
        (atropos, pad_segm, [("classified_image", "op1")]),
        (pad_segm, sel_labels, [("output_image", "in_segm")]),
        (sel_labels, get_wm, [("out_wm", "op1")]),

def _list_outputs(self):
        outputs = super(ImageMath, self)._list_outputs()
        if self.inputs.copy_header:  # Fix headers
            _copy_header(self.inputs.op1, outputs["output_image"], set_dtype=False)
        return outputs

)
    me_csf = pe.Node(ImageMath(operation="ME", op2="10"), name="10_me_csf")

    # ImageMath ${DIMENSION} ${EXTRACTION_GM} addtozero ${EXTRACTION_GM} ${EXTRACTION_TMP}
    # MultiplyImages ${DIMENSION} ${EXTRACTION_GM} ${ATROPOS_GM_CLASS_LABEL} ${EXTRACTION_GM}
    # ImageMath ${DIMENSION} ${EXTRACTION_SEGMENTATION} addtozero ${EXTRACTION_WM} ${EXTRACTION_GM}
    add_gm = pe.Node(ImageMath(operation="addtozero"), name="11_add_gm")
    relabel_gm = pe.Node(
        MultiplyImages(
            dimension=3,
            second_input=in_segmentation_model[-2],
            output_product_image="12_relabel_gm.nii.gz",
        ),
        name="12_relabel_gm",
    )
    add_gm_wm = pe.Node(ImageMath(operation="addtozero"), name="13_add_gm_wm")

    # Superstep 7
    # Split segmentation in binary masks
    sel_labels2 = pe.Node(
        niu.Function(function=_select_labels, output_names=["out_gm", "out_wm"]),
        name="14_sel_labels2",
    )
    sel_labels2.inputs.labels = in_segmentation_model[2:]

    # ImageMath ${DIMENSION} ${EXTRACTION_MASK} addtozero ${EXTRACTION_MASK} ${EXTRACTION_TMP}
    add_7 = pe.Node(ImageMath(operation="addtozero"), name="15_add_7")
    # ImageMath ${DIMENSION} ${EXTRACTION_MASK} ME ${EXTRACTION_MASK} 2
    me_7 = pe.Node(ImageMath(operation="ME", op2="2"), name="16_me_7")
    # ImageMath ${DIMENSION} ${EXTRACTION_MASK} GetLargestComponent ${EXTRACTION_MASK}
    comp_7 = pe.Node(ImageMath(operation="GetLargestComponent"), name="17_comp_7")
    # ImageMath ${DIMENSION} ${EXTRACTION_MASK} MD ${EXTRACTION_MASK} 4

number_of_tissue_classes=in_segmentation_model[0],
            n_iterations=3,
            convergence_threshold=0.0,
            mrf_radius=[1, 1, 1],
            mrf_smoothing_factor=0.1,
            likelihood_model="Gaussian",
            use_random_seed=use_random_seed,
        ),
        name="01_atropos",
        n_procs=omp_nthreads,
        mem_gb=mem_gb,
    )

    # massage outputs
    pad_segm = pe.Node(
        ImageMath(operation="PadImage", op2="%d" % padding), name="02_pad_segm"
    )
    pad_mask = pe.Node(
        ImageMath(operation="PadImage", op2="%d" % padding), name="03_pad_mask"
    )

    # Split segmentation in binary masks
    sel_labels = pe.Node(
        niu.Function(
            function=_select_labels, output_names=["out_wm", "out_gm", "out_csf"]
        ),
        name="04_sel_labels",
    )
    sel_labels.inputs.labels = list(reversed(in_segmentation_model[1:]))

    # Select largest components (GM, WM)
    # ImageMath ${DIMENSION} ${EXTRACTION_WM} GetLargestComponent ${EXTRACTION_WM}

"bias_image",
                "out_segm",
                "out_tpms",
            ]
        ),
        name="outputnode",
    )

    copy_xform = pe.Node(
        CopyXForm(fields=["out_file", "out_mask", "bias_corrected", "bias_image"]),
        name="copy_xform",
        run_without_submitting=True,
    )

    trunc = pe.MapNode(
        ImageMath(operation="TruncateImageIntensity", op2="0.01 0.999 256"),
        name="truncate_images",
        iterfield=["op1"],
    )
    inu_n4 = pe.MapNode(
        N4BiasFieldCorrection(
            dimension=3,
            save_bias=False,
            copy_header=True,
            n_iterations=[50] * 4,
            convergence_threshold=1e-7,
            shrink_factor=4,
            bspline_fitting_distance=bspline_fitting_distance,
        ),
        n_procs=omp_nthreads,
        name="inu_n4",
        iterfield=["input_image"],

md_7_2 = pe.Node(ImageMath(operation="MD", op2="5"), name="21_md_7_2")
    # ImageMath ${DIMENSION} ${EXTRACTION_MASK} ME ${EXTRACTION_MASK} 5
    me_7_2 = pe.Node(ImageMath(operation="ME", op2="5"), name="22_me_7_2")

    # De-pad
    depad_mask = pe.Node(
        ImageMath(operation="PadImage", op2="-%d" % padding), name="23_depad_mask"
    )
    depad_segm = pe.Node(
        ImageMath(operation="PadImage", op2="-%d" % padding), name="24_depad_segm"
    )
    depad_gm = pe.Node(
        ImageMath(operation="PadImage", op2="-%d" % padding), name="25_depad_gm"
    )
    depad_wm = pe.Node(
        ImageMath(operation="PadImage", op2="-%d" % padding), name="26_depad_wm"
    )
    depad_csf = pe.Node(
        ImageMath(operation="PadImage", op2="-%d" % padding), name="27_depad_csf"
    )

    msk_conform = pe.Node(niu.Function(function=_conform_mask), name="msk_conform")
    merge_tpms = pe.Node(niu.Merge(in_segmentation_model[0]), name="merge_tpms")
    # fmt: off
    wf.connect([
        (inputnode, copy_xform, [(("in_files", _pop), "hdr_file")]),
        (inputnode, pad_mask, [("in_mask", "op1")]),
        (inputnode, atropos, [
            ("in_files", "intensity_images"),
            ("in_mask_dilated", "mask_image"),
        ]),
        (inputnode, msk_conform, [(("in_files", _pop), "in_reference")]),

output_product_image="12_relabel_gm.nii.gz",
        ),
        name="12_relabel_gm",
    )
    add_gm_wm = pe.Node(ImageMath(operation="addtozero"), name="13_add_gm_wm")

    # Superstep 7
    # Split segmentation in binary masks
    sel_labels2 = pe.Node(
        niu.Function(function=_select_labels, output_names=["out_gm", "out_wm"]),
        name="14_sel_labels2",
    )
    sel_labels2.inputs.labels = in_segmentation_model[2:]

    # ImageMath ${DIMENSION} ${EXTRACTION_MASK} addtozero ${EXTRACTION_MASK} ${EXTRACTION_TMP}
    add_7 = pe.Node(ImageMath(operation="addtozero"), name="15_add_7")
    # ImageMath ${DIMENSION} ${EXTRACTION_MASK} ME ${EXTRACTION_MASK} 2
    me_7 = pe.Node(ImageMath(operation="ME", op2="2"), name="16_me_7")
    # ImageMath ${DIMENSION} ${EXTRACTION_MASK} GetLargestComponent ${EXTRACTION_MASK}
    comp_7 = pe.Node(ImageMath(operation="GetLargestComponent"), name="17_comp_7")
    # ImageMath ${DIMENSION} ${EXTRACTION_MASK} MD ${EXTRACTION_MASK} 4
    md_7 = pe.Node(ImageMath(operation="MD", op2="4"), name="18_md_7")
    # ImageMath ${DIMENSION} ${EXTRACTION_MASK} FillHoles ${EXTRACTION_MASK} 2
    fill_7 = pe.Node(ImageMath(operation="FillHoles", op2="2"), name="19_fill_7")
    # ImageMath ${DIMENSION} ${EXTRACTION_MASK} addtozero ${EXTRACTION_MASK} \
    # ${EXTRACTION_MASK_PRIOR_WARPED}
    add_7_2 = pe.Node(ImageMath(operation="addtozero"), name="20_add_7_2")
    # ImageMath ${DIMENSION} ${EXTRACTION_MASK} MD ${EXTRACTION_MASK} 5
    md_7_2 = pe.Node(ImageMath(operation="MD", op2="5"), name="21_md_7_2")
    # ImageMath ${DIMENSION} ${EXTRACTION_MASK} ME ${EXTRACTION_MASK} 5
    me_7_2 = pe.Node(ImageMath(operation="ME", op2="5"), name="22_me_7_2")