Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

[MDAKit/ Hackathon] Pocket Volume calculation #48

Open
H-EKE opened this issue Sep 29, 2023 · 1 comment
Open

[MDAKit/ Hackathon] Pocket Volume calculation #48

H-EKE opened this issue Sep 29, 2023 · 1 comment

Comments

@H-EKE
Copy link

H-EKE commented Sep 29, 2023

Hi,

We are trying to implement Fpocket (https://github.com/Discngine/fpocket) into MDAnalysis

@H-EKE
Copy link
Author

H-EKE commented Sep 29, 2023

###DRAFT

We are starting by adapting the same procedure as with HOLE for Fpocket :)



def run_fpocket(outfile, infile_text, executable):
    """Run the HOLE program.

    Parameters
    ----------
    outfile: str
        Output file name
    infile_text: str
        HOLE input text
        (typically generated by :func:`set_up_hole_input`)
    executable: str
        HOLE executable

    Returns
    -------
    str
        Output file name
    """
    with open(outfile, 'w') as output:
        proc = subprocess.Popen(executable, stdin=subprocess.PIPE,
                                stdout=output)
        stdout, stderr = proc.communicate(infile_text.encode('utf-8'))

    # check output in case of errors
    with open(outfile, 'r') as output:
        for line in output:
            if line.strip().startswith(('*** ERROR ***', 'ERROR')):
                proc.returncode = 255
                break

    # die in case of error
    if proc.returncode != 0:
        msg = 'HOLE failure ({}). Check output {}'
        logger.fatal(msg.format(proc.returncode, outfile))
        if stderr is not None:
            logger.fatal(stderr)
        raise ApplicationError(proc.returncode,
                               msg.format(executable, outfile))

    logger.info('HOLE finished. Output: {}'.format(outfile))
    return outfile
    def fpocket(pdbfile,
         outfile='hole.out',
         executable='hole',
         tmpdir=os.path.curdir,
         sample=0.2,
         end_radius=22.0,
         cpoint=None,
         cvect=None,
         random_seed=None,
         ignore_residues=IGNORE_RESIDUES,
         output_level=0,
         dcd=None,
         dcd_iniskip=0,
         dcd_step=1,
         keep_files=True):
         
def collect_fpocket(outfile='hole.out'):
    """Collect data from HOLE output.

    Parameters
    ----------
    outfile: str, optional
        HOLE output file to read. Default: 'hole.out'


    Returns
    -------
    dict
        Dictionary of HOLE profiles as record arrays
    """
    fmt = util.FORTRANReader('3F12')
    recarrays = {}

    with open(outfile, 'r') as output:
        toggle_read = False
        profile = 0
        records = []
        for line in output:
            line = line.rstrip()  # preserve columns in FORTRAN output
            # multiple frames from dcd in?
            if line.startswith(" Starting calculation for position number"):
                fields = line.split()
                profile = int(fields[5])
                records = []
                logger.debug('Started reading profile {}'.format(profile))
                continue

            # found data
            if line.startswith(' cenxyz.cvec'):
                toggle_read = True
                logger.debug('Started reading data')
                continue

            if toggle_read:
                if len(line.strip()) != 0:
                    try:
                        rxncoord, radius, cenlineD = fmt.read(line)
                    except:
                        msg = 'Problem parsing line: {}. Check output file {}'
                        logger.exception(msg.format(line, outfile))
                        raise
                    records.append((rxncoord, radius, cenlineD))
                    continue
                # end of data
                else:
                    toggle_read = False
                    names = ['rxn_coord', 'radius', 'cen_line_D']
                    recarr = np.rec.fromrecords(records,
                                                formats='f8, f8, f8',
                                                names=names)
                    recarrays[profile] = recarr

        return recarrays

profiles = hole2.hole(PDB_HOLE, executable='~/hole2/exe/hole')

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

1 participant