How to use the parmed.amber.AmberParm function in ParmEd

def __init__(self, constraints=app.HBonds, rigid_water=True, nonbondedCutoff=8.0 * unit.angstroms, use_dispersion_correction=True, nonbondedMethod=app.PME, hydrogenMass=None, switch_width=None, ewaldErrorTolerance=5E-4, **kwargs):

        TestSystem.__init__(self, **kwargs)

        try:
            from parmed.amber import AmberParm
        except ImportError as e:
            print("DHFR test system requires Parmed (`import parmed`).")
            raise(e)

        prmtop_filename = get_data_filename("data/dhfr/prmtop")
        crd_filename = get_data_filename("data/dhfr/inpcrd")

        # Initialize system.
        self.prmtop = AmberParm(prmtop_filename, crd_filename)
        system = self.prmtop.createSystem(constraints=constraints, nonbondedMethod=nonbondedMethod, rigidWater=rigid_water, nonbondedCutoff=nonbondedCutoff, hydrogenMass=hydrogenMass)

        # Extract topology
        self.topology = self.prmtop.topology

        # Set dispersion correction use.
        forces = {system.getForce(index).__class__.__name__: system.getForce(index) for index in range(system.getNumForces())}
        forces['NonbondedForce'].setUseDispersionCorrection(use_dispersion_correction)
        forces['NonbondedForce'].setEwaldErrorTolerance(ewaldErrorTolerance)

        if switch_width is not None:
            forces['NonbondedForce'].setUseSwitchingFunction(True)
            forces['NonbondedForce'].setSwitchingDistance(nonbondedCutoff - switch_width)


        positions = self.prmtop.positions

assert isinstance(precision, int), "Precision %s is not an integer." % precision

    #Import ParmEd
    import parmed
    #Require version 2.0.4 or later of ParmEd, otherwise ParmEd corrupts [ defaults ] section in GROMACS topologies with incorrect FudgeLJ/FudgeQQ
    try:
        ver = parmed.version
    except:
        oldParmEd = Exception('ERROR: ParmEd is too old, please upgrade to 2.0.4 or later')
        raise oldParmEd
    if ver < (2,0,4):
        raise RuntimeError("ParmEd is too old, please upgrade to 2.0.4 or later")


    #Read AMBER to ParmEd object
    structure = parmed.amber.AmberParm( in_prmtop, in_crd )
    #Make GROMACS topology
    gromacs_topology = parmed.gromacs.GromacsTopologyFile.from_structure( structure )
    #Write
    parmed.gromacs.GromacsTopologyFile.write( gromacs_topology, out_top )
    if precision == None:
        parmed.gromacs.GromacsGroFile.write( gromacs_topology, out_gro )
    else:
        parmed.gromacs.GromacsGroFile.write( gromacs_topology, out_gro, precision = precision )

    return out_top, out_gro

# Filter titratable residues based on min and max pKa
    new_reslist = []
    for res in titratable_residues:
        if getattr(residues, res).Eo < mineo: continue
        if getattr(residues, res).Eo > maxeo: continue
        new_reslist.append(res)
    titratable_residues = new_reslist
    del new_reslist

    # Make sure we still have a couple residues
    if len(titratable_residues) == 0:
        raise AmberError('No titratable residues fit your criteria!')

    # Load the topology file
    parm = AmberParm(opt.prmtop)

    # Replace an un-set notresnums with an empty list so we get __contains__()
    if notresnums is None:
        notresnums = []

    # If we have a list of residue numbers, check that they're all titratable
    if resnums is not None:
        for resnum in resnums:
            if resnum > parm.ptr('nres'):
                raise AmberError('%s only has %d residues. (%d chosen)' %
                                     (parm, parm.ptr('nres'), resnum))
            if resnum <= 0:
                raise AmberError('Cannot select negative residue numbers.')
            resname = parm.parm_data['RESIDUE_LABEL'][resnum-1]
            if not resname in titratable_residues:
                raise AmberError('Residue number %s [%s] is not titratable'

if not resname in residues.titratable_residues:
                raise AmberError('%s is not a titratable residue!' % resname)
            titratable_residues.append(resname)
    else:
        for resname in residues.titratable_residues:
            titratable_residues.append(resname)

    solvent_ions = ['WAT', 'Na+', 'Br-', 'Cl-', 'Cs+', 'F-', 'I-', 'K+', 'Li+',
                    'Mg+', 'Rb+', 'CIO', 'IB', 'MG2']

    # Make sure we still have a couple residues
    if len(titratable_residues) == 0:
        raise AmberError('No titratable residues fit your criteria!')

    # Load the topology file
    parm = AmberParm(opt.prmtop)

    # Replace an un-set notresnums with an empty list so we get __contains__()
    if notresnums is None:
        notresnums = []

    # If we have a list of residue numbers, check that they're all titratable
    if resnums is not None:
        for resnum in resnums:
            if resnum > parm.ptr('nres'):
                raise AmberError('%s only has %d residues. (%d chosen)' %
                                     (parm, parm.ptr('nres'), resnum))
            if resnum <= 0:
                raise AmberError('Cannot select negative residue numbers.')
            resname = parm.parm_data['RESIDUE_LABEL'][resnum-1]
            if not resname in titratable_residues:
                raise AmberError('Residue number %s [%s] is not titratable'

def default(self, obj):
        """If input object is an ndarray it will be converted into a dict
        holding dtype, shape and the data, base64 encoded.
        """
        if isinstance(obj, AmberParm):
            log.info("Encountered AmberParm, returning name.")
            return obj.name
        if isinstance(obj, Structure):
            log.warning("Encountered Structure, which does not store filename.")
            return ""

        if isinstance(obj, np.ndarray):
            if obj.flags["C_CONTIGUOUS"]:
                obj_data = obj.data
            else:
                cont_obj = np.ascontiguousarray(obj)
                assert cont_obj.flags["C_CONTIGUOUS"]
                obj_data = cont_obj.data
            data_b64 = base64.b64encode(obj_data)
            # obj_data = obj.tolist()
            return dict(

# Filter titratable residues based on min and max pKa
    new_reslist = []
    for res in titratable_residues:
        if getattr(residues, res).pKa < minpka: continue
        if getattr(residues, res).pKa > maxpka: continue
        new_reslist.append(res)
    titratable_residues = new_reslist
    del new_reslist

    # Make sure we still have a couple residues
    if len(titratable_residues) == 0:
        raise AmberError('No titratable residues fit your criteria!')

    # Load the topology file
    parm = AmberParm(opt.prmtop)

    # Replace an un-set notresnums with an empty list so we get __contains__()
    if notresnums is None:
        notresnums = []

    # If we have a list of residue numbers, check that they're all titratable
    if resnums is not None:
        for resnum in resnums:
            if resnum > parm.ptr('nres'):
                raise AmberError('%s only has %d residues. (%d chosen)' %
                                     (parm, parm.ptr('nres'), resnum))
            if resnum <= 0:
                raise AmberError('Cannot select negative residue numbers.')
            resname = parm.parm_data['RESIDUE_LABEL'][resnum-1]
            if not resname in titratable_residues:
                raise AmberError('Residue number %s [%s] is not titratable'

def get_rmsd(initparas):

    global idxs, mcresids2, atompairs

    #Modify the C4 terms in the prmtop file
    for i in range(0, len(idxs)):
        prmtop.parm_data['LENNARD_JONES_CCOEF'][idxs[i]] = initparas[i]

    #Overwrite the prmtop file
    prmtop.write_parm('OptC4.top')

    #Perform the OpenMM optimization
    #Use AmberParm function to transfer the topology and
    #coordinate file to the object OpenMM can use
    Ambermol = AmberParm('OptC4.top', options.cfile)

    # Create the OpenMM system
    print('Creating OpenMM System')
    if options.simupha == 'gas':
        system = Ambermol.createSystem(nonbondedMethod=app.NoCutoff)
    elif options.simupha == 'liquid':
        system = Ambermol.createSystem(nonbondedMethod=app.PME,
                                       nonbondedCutoff=8.0*u.angstroms,
                                       constraints=app.HBonds,)

    #Add restraints
    force = mm.CustomExternalForce("k*((x-x0)^2+(y-y0)^2+(z-z0)^2)")
    force.addGlobalParameter("k", 200.0)
    force.addPerParticleParameter("x0")
    force.addPerParticleParameter("y0")
    force.addPerParticleParameter("z0")