Python stages.Stages类代码示例

def nlin_part(xfm : XfmHandler, inv_xfm : Optional[XfmHandler] = None) -> Result[XfmHandler]:
    """
    *** = non linear deformations
    --- = linear (affine) deformations

    Input:
    xfm     :     ******------>
    inv_xfm :    <******------ (optional)

    Calculated:
    inv_lin_xfm :      <------

    Returned:
    concat :      ******------> +
                       <------
    equals :      ******>

    Compute the nonlinear part of a transform as follows:
    go forwards across xfm and then backwards across the linear part
    of the inverse xfm (by first calculating the inverse or using the one supplied) 
    Finally, use minc_displacement to compute the resulting gridfile of the purely 
    nonlinear part.

    The optional inv_xfm (which must be the inverse!) is an optimization -
    we don't go looking for an inverse by filename munging and don't programmatically
    keep a log of operations applied, so any preexisting inverse must be supplied explicitly.
    """
    s = Stages()
    inv_xfm = inv_xfm or s.defer(invert_xfmhandler(xfm))
    inv_lin_part = s.defer(lin_from_nlin(inv_xfm)) 
    xfm = s.defer(concat_xfmhandlers([xfm, inv_lin_part]))
    return Result(stages=s, output=xfm)

def mincblob(op : str, grid : MincAtom, subdir : str = "tmp") -> Result[MincAtom]:
    """
    Low-level mincblob wrapper with the one exception being the determinant option. By
    default the inner clockwork of mincblob subtracts 1 from all determinant values that
    are being calculated. As such, 1 needs to be added to the result of the mincblob call.
    We will do that here, because it makes most sense here.
    >>> stages = mincblob('determinant', MincAtom("/images/img_1.mnc", pipeline_sub_dir="/tmp")).stages
    >>> [s.render() for s in stages]
    ['mincblob -clobber -determinant /images/img_1.mnc /tmp/img_1/img_1_determinant.mnc']
    """
    if op not in ["determinant", "trace", "translation", "magnitude"]:
        raise ValueError('mincblob: invalid operation %s' % op)

    # if we are calculating the determinant, the first file produced is a temp file:
    if op == "determinant":
        out_file = grid.newname_with_suffix("_temp_det", subdir=subdir)
    else:
        out_file = grid.newname_with_suffix('_' + op, subdir=subdir)

    stage = CmdStage(inputs=(grid,), outputs=(out_file,),
                 cmd=['mincblob', '-clobber', '-' + op, grid.path, out_file.path])

    s = Stages([stage])
    # now create the proper determinant if that's what was asked for
    if op == "determinant":
        result_file = s.defer(mincmath(op='add',
                                       const=1,
                                       vols=[out_file],
                                       subdir=subdir,
                                       new_name=grid.filename_wo_ext + "_det"))
    else:
        result_file = out_file

    return Result(stages=s, output=result_file)

def lsq6_pipeline(options):
    # TODO could also allow pluggable pipeline parts e.g. LSQ6 could be substituted out for the modified LSQ6
    # for the kidney tips, etc...
    output_dir    = options.application.output_directory
    pipeline_name = options.application.pipeline_name

    # TODO this is tedious and annoyingly similar to the registration chain and MBM ...
    lsq6_dir      = os.path.join(output_dir, pipeline_name + "_lsq6")
    processed_dir = os.path.join(output_dir, pipeline_name + "_processed")
    imgs = get_imgs(options.application)

    s = Stages()

    # TODO this is quite tedious and duplicates stuff in the registration chain ...
    resolution = (options.registration.resolution or
                  get_resolution_from_file(
                      s.defer(registration_targets(lsq6_conf=options.lsq6,
                                                   app_conf=options.application,
                                                   reg_conf=options.registration)).registration_standard.path))

    # FIXME: why do we have to call registration_targets *outside* of lsq6_nuc_inorm? is it just because of the extra
    # options required?
    targets = s.defer(registration_targets(lsq6_conf=options.lsq6,
                                   app_conf=options.application,
                                   reg_conf=options.registration,
                                   first_input_file=imgs[0]))
    # This must happen after calling registration_targets otherwise it will resample to options.registration.resolution
    options.registration = options.registration.replace(resolution=resolution)
    lsq6_result = s.defer(lsq6_nuc_inorm(imgs=imgs,
                                         resolution=resolution,
                                         registration_targets=targets,
                                         lsq6_dir=lsq6_dir,
                                         lsq6_options=options.lsq6))

    return Result(stages=s, output=lsq6_result)

def asymmetry_pipeline(options):

    output_dir    = options.application.output_directory
    pipeline_name = options.application.pipeline_name
    processed_dir = os.path.join(output_dir, pipeline_name + "_processed")

    s = Stages()

    #imgs_ = [MincAtom(f, pipeline_sub_dir=processed_dir) for f in options.application.files]

    imgs_ = get_imgs(options.application)

    check_MINC_input_files([img.path for img in imgs_])

    imgs  = pd.Series(imgs_, index=[img.filename_wo_ext for img in imgs_])
    flipped_imgs = imgs.apply(lambda img: s.defer(volflip(img)))  # TODO add flags to control flip axis ...

    # TODO ugly - MincAtom API should allow this somehow without mutation (also, how to pass into `volflip`, etc.?)
    for f_i in flipped_imgs:
        f_i.output_sub_dir += "_flipped"

    check_MINC_input_files(imgs.apply(lambda img: img.path))

    grouped_files_df = pd.DataFrame({'file' : pd.concat([imgs, flipped_imgs])}).assign(group=lambda df: df.index)

    two_level_result = s.defer(two_level(grouped_files_df, options=options))

    return Result(stages=s, output=two_level_result)

def cortical_thickness_pipeline(options):
    s = Stages()

    #imgs = [MincAtom(name, pipeline_sub_dir=os.path.join(options.application.output_directory,
    #                                                     options.application.pipeline_name + "_processed"))
    #        for name in options.application.files]

    pipeline_sub_dir = os.path.join(options.application.output_directory,
                                    options.application.pipeline_name + "_processed")

    #def atom(atom_type, file):
    #    return atom_type(file, pipeline_sub_dir=pipeline_sub_dir)  # TODO output_sub_dir, ....

    # TODO are all these fields actually used?  If not, omit from CSV?
    xfms = (pd.read_csv(options.thickness.xfm_csv)
            .apply(axis=1,  # TODO fill out <..>Atom(...) fields ...
                   func=lambda row: XfmHandler(
                          source=MincAtom(row.source, pipeline_sub_dir=pipeline_sub_dir),
                          target=MincAtom(row.target, pipeline_sub_dir=pipeline_sub_dir),
                          resampled=None,   #MincAtom(row.resampled, pipeline_sub_dir=pipeline_sub_dir),
                          xfm=XfmAtom(row.xfm, pipeline_sub_dir=pipeline_sub_dir))))
    # TODO better way to unpack?

    result = s.defer(cortical_thickness(xfms=xfms,
                                        atlas=NotImplemented,
                                        label_mapping=options.thickness.label_mapping,
                                        atlas_fwhm=options.thickness.atlas_fwhm,
                                        thickness_fwhm=options.thickness.thickness_fwhm))

    return Result(stages=s, output=result)

def tamarack_pipeline(options):

    output_dir    = options.application.output_directory
    pipeline_name = options.application.pipeline_name
    #processed_dir = os.path.join(output_dir, pipeline_name + "_processed")
    first_level_dir = os.path.join(output_dir, pipeline_name + "_first_level")

    s = Stages()

    with open(options.application.csv_file, 'r') as f:
        files_df = (pd.read_csv(filepath_or_buffer=f,
                                usecols=['group', 'filename'])
                    .assign(file=lambda df:
                                   df.apply(axis="columns",
                                            func=lambda r:
                                                   MincAtom(r.filename.strip(),
                                                            pipeline_sub_dir=os.path.join(first_level_dir,
                                                                                          "%s_processed" % r.group.strip())))))

    check_MINC_input_files(files_df.file.apply(lambda img: img.path))

    #grouped_files_df = pd.DataFrame({'file' : pd.concat([imgs])}).assign(group=lambda df: df.index)

    tamarack_result = s.defer(tamarack(files_df, options=options))

    tamarack_result.first_level_results.applymap(maybe_deref_path).to_csv("first_level_results.csv", index=False)
    tamarack_result.resampled_determinants.applymap(maybe_deref_path).to_csv("resampled_determinants.csv", index=False)
    tamarack_result.overall_determinants.applymap(maybe_deref_path).to_csv("overall_determinants.csv", index=False)

    return Result(stages=s, output=tamarack_result)

 def scale_transform(xfm, scale, newname_wo_ext):
     s = Stages()
     defs = s.defer(as_deformation(transform=xfm.xfm, reference=xfm.source))
     scaled_defs = (defs.xfm.newname(newname_wo_ext) if newname_wo_ext else
                     defs.xfm.newname_with_suffix("_scaled_%s" % scale))
     s.defer(CmdStage(cmd=['c3d', '-scale', str(scale), defs.path, "-o", scaled_defs.path],
                      inputs=(defs,), outputs=(scaled_defs,)))
     return Result(stages=s, output=scaled_defs)

def determinant(displacement_grid : MincAtom) -> Result[MincAtom]:
    """
    Takes a displacement field (deformation grid, vector field, those are
    all the same thing) and calculates the proper determinant (mincblob()
    takes care of adding 1 to the silly output of running mincblob directly)
    """
    s = Stages()
    det = s.defer(mincblob(op='determinant', grid=displacement_grid))
    return Result(stages=s, output=det)

def convert(infile : ImgAtom, out_ext : str) -> Result[ImgAtom]:
    s = Stages()
    outfile = infile.newext(ext=out_ext)
    if infile.mask is not None:
        outfile.mask = s.defer(convert(infile.mask, out_ext=out_ext))
    if infile.labels is not None:
        outfile.mask = s.defer(convert(infile.labels, out_ext=out_ext))
    s.add(CmdStage(inputs=(infile,), outputs=(outfile,),
                   cmd = ['c3d', infile.path, '-o', outfile.path]))
    return Result(stages=s, output=outfile)

 def to_mni_xfm(xfm):
     s = Stages()
     defs = xfm.newname_with_suffix("_defs", subdir="tmp")
     s.add(CmdStage(cmd=["transformix", "-def", "all",
                         "-out", defs.dir,
                         "-tp", xfm.path,
                         "-xfm", os.path.join(defs.filename_wo_ext, defs.ext)],
                    inputs=(xfm,), outputs=(defs,)))
     out_xfm = s.defer(itk.itk_convert_xfm(defs, out_ext=".mnc"))
     return Result(stages=s, output=out_xfm)

 def average_transforms(xfms, avg_xfm):
     s = Stages()
     defs = [s.defer(as_deformation(transform=xfm.xfm, reference_image=xfm.source)) for xfm in xfms]
     #avg_img = NotImplemented
     avg = imageToXfm(s.defer(average_images(defs,
                                             avg_file=xfmToImage(avg_xfm),
                                             #output_dir=os.path.join(defs[0].pipeline_sub_dir,
                                             #                        defs[0].output_sub_dir,
                                             #                        "transforms")
                                             )))
     return Result(stages=s, output=avg)

def NLIN_pipeline(options):

    # if options.application.files is None:
    #     raise ValueError("Please, some files! (or try '--help')")  # TODO make a util procedure for this

    output_dir    = options.application.output_directory
    pipeline_name = options.application.pipeline_name

    # TODO this is tedious and annoyingly similar to the registration chain and MBM and LSQ6 ...
    processed_dir = os.path.join(output_dir, pipeline_name + "_processed")
    nlin_dir      = os.path.join(output_dir, pipeline_name + "_nlin")

    resolution = (options.registration.resolution  # TODO does using the finest resolution here make sense?
                  or min([get_resolution_from_file(f) for f in options.application.files]))

    imgs = get_imgs(options.application)

    # imgs = [MincAtom(f, pipeline_sub_dir=processed_dir) for f in options.application.files]

    # determine NLIN settings by overriding defaults with
    # any settings present in protocol file, if it exists
    # could add a hook to print a message announcing completion, output files,
    # add more stages here to make a CSV

    initial_target_mask = MincAtom(options.nlin.target_mask) if options.nlin.target_mask else None
    initial_target = MincAtom(options.nlin.target, mask=initial_target_mask)

    full_hierarchy = get_nonlinear_configuration_from_options(nlin_protocol=options.nlin.nlin_protocol,
                                                              reg_method=options.nlin.reg_method,
                                                              file_resolution=resolution)

    s = Stages()

    nlin_result = s.defer(nlin_build_model(imgs, initial_target=initial_target, conf=full_hierarchy, nlin_dir=nlin_dir))

    # TODO return these?
    inverted_xfms = [s.defer(invert_xfmhandler(xfm)) for xfm in nlin_result.output]

    if options.stats.calc_stats:
        # TODO: put the stats part behind a flag ...

        determinants = [s.defer(determinants_at_fwhms(
                                  xfm=inv_xfm,
                                  inv_xfm=xfm,
                                  blur_fwhms=options.stats.stats_kernels))
                        for xfm, inv_xfm in zip(nlin_result.output, inverted_xfms)]

        return Result(stages=s,
                      output=Namespace(nlin_xfms=nlin_result,
                                       avg_img=nlin_result.avg_img,
                                       determinants=determinants))
    else:
        # there's no consistency in what gets returned, yikes ...
        return Result(stages=s, output=Namespace(nlin_xfms=nlin_result, avg_img=nlin_result.avg_img))

 def f(imgs: List[MincAtom],
       initial_target: MincAtom,
       conf: reg_module.Conf,
       nlin_dir: str,
       nlin_prefix: str,
       tournament_name_wo_ext: str = "tournament") -> Result[List[XfmHandler]]:
   s = Stages()
   Weight = int
   def h(xfms : List[XfmHandler], name_wo_ext : str) -> List[XfmHandler]:
       # TODO add weights to each return
       # TODO check len(...) == 0 case??
       if len(xfms) <= 1:
           return xfms
       else:
           first_half  = xfms[: len(xfms)//2]
           second_half = xfms[len(xfms)//2 :]
           first_half_result  = h(first_half,  name_wo_ext=name_wo_ext + "_L")
           second_half_result = h(second_half, name_wo_ext=name_wo_ext + "_R")
           A_halfway_to_B, B_halfway_to_A, avg_img = s.defer(
               nonlinear_midpoint_xfm(
                   img_A = first_half_result[0].resampled,
                   img_B = second_half_result[0].resampled,
                   out_name_wo_ext=name_wo_ext,
                   nlin_algorithm=reg_module,
                   conf=conf,
                   out_dir=nlin_dir))
           xfms_to_midpoint = ([XfmHandler(source=xfm.source,
                                           target=avg_img,
                                           resampled=A_halfway_to_B.resampled,
                                           xfm=s.defer(xfmconcat([xfm.xfm, A_halfway_to_B.xfm],
                                                                 name="%s_%s" % (xfm.source.filename_wo_ext,
                                                                                 name_wo_ext))))
                                for xfm in first_half_result]
                               + [XfmHandler(source=xfm.source,
                                             target=avg_img,
                                             resampled=B_halfway_to_A.resampled,
                                             xfm=s.defer(xfmconcat([xfm.xfm, B_halfway_to_A.xfm],
                                                                   name="%s_%s" % (xfm.source.filename_wo_ext,
                                                                                   name_wo_ext))))
                                  for xfm in second_half_result])
           return xfms_to_midpoint
   identity_xfm = s.defer(param2xfm(out_xfm=XfmAtom(pipeline_sub_dir=imgs[0].pipeline_sub_dir,
                                                    output_sub_dir=imgs[0].output_sub_dir,
                                                    name=os.path.join(imgs[0].pipeline_sub_dir,
                                                                      imgs[0].output_sub_dir,
                                                                      "id.xfm"))))
   initial_xfms = [XfmHandler(source=img, target=img,
                              resampled=img, xfm=identity_xfm) for img in imgs]
   xfms_to_avg = h(initial_xfms, tournament_name_wo_ext)
   avg_img = xfms_to_avg[0].target
   return Result(stages=s, output=WithAvgImgs(avg_img=avg_img, avg_imgs=[avg_img],
                                              output=xfms_to_avg))

 def average_transforms(xfms, avg_xfm):
     intermediate_xfm = avg_xfm.newname_with_suffix("_inter", subdir="tmp")
     s = Stages()
     s.add(CmdStage(cmd=["echo", ('(Transform "WeightedCombinationTransform")\n'
                                  '(SubTransforms %s)\n'
                                  '(NormalizeCombinationsWeights "true")\n') %
                                    ' '.join(sorted(xfm.path for xfm in xfms))],
                    inputs=xfms, outputs=(intermediate_xfm,)))
     s.add(CmdStage(cmd=["transformix", "-def", "all",
                         "-out", os.path.dirname(avg_xfm.path),
                         "-tp", intermediate_xfm.path,
                         "-xfm", avg_xfm.path],
                    inputs=(intermediate_xfm,), outputs=(avg_xfm,)))

def nlin_displacement(xfm : XfmHandler, inv_xfm : Optional[XfmHandler] = None) -> Result[MincAtom]:
    """
    See: nlin_part().

    This returns the nonlinear part of the input
    transformation (xfm) in the form of a grid file (vector field).
    All transformations are encapsulated in this field (linear parts
    that are normally specified in the .xfm file are placed in the
    vector field)
    """
    
    s = Stages()
    return Result(stages=s,
                  output=s.defer(minc_displacement(
                                   s.defer(nlin_part(xfm, inv_xfm=inv_xfm)))))