text:
/*************************************************************************
 * GLFW 3.2 - www.glfw.org
 * A library for OpenGL, window and input
 *------------------------------------------------------------------------
 * Copyright (c) 2002-2006 Marcus Geelnard
 * Copyright (c) 2006-2016 Camilla Berglund <elmindreda@glfw.org>
 *
 * This software is provided 'as-is', without any express or implied
 * warranty. In no event will the authors be held liable for any damages
 * arising from the use of this software.
 *
 * Permission is granted to anyone to use this software for any purpose,
 * including commercial applications, and to alter it and redistribute it
 * freely, subject to the following restrictions:
 *
 * 1. The origin of this software must not be misrepresented; you must not
 *    claim that you wrote the original software. If you use this software
 *    in a product, an acknowledgment in the product documentation would
 *    be appreciated but is not required.
 *
 * 2. Altered source versions must be plainly marked as such, and must not
 *    be misrepresented as being the original software.
 *
 * 3. This notice may not be removed or altered from any source
 *    distribution.
 *
 *************************************************************************/

#ifndef _glfw3_native_h_
#define _glfw3_native_h_

#ifdef __cplusplus
extern "C" {
#endif

/*************************************************************************
 * Doxygen documentation
 *************************************************************************/

/*! @file glfw3native.h
 *  @brief The header of the native access functions.
 *
 *  This is the header file of the native access functions.  See @ref native for
 *  more information.
 */
/*! @defgroup native Native access
 *
 *  **By using the native access functions you assert that you know what you're
 *  doing and how to fix problems caused by using them.  If you don't, you
 *  shouldn't be using them.**
 *
 *  Before the inclusion of @ref glfw3native.h, you may define exactly one
 *  window system API macro and zero or more context creation API macros.
 *
 *  The chosen backends must match those the library was compiled for.  Failure
 *  to do this will cause a link-time error.
 *
 *  The available window API macros are:
 *  * `GLFW_EXPOSE_NATIVE_WIN32`
 *  * `GLFW_EXPOSE_NATIVE_COCOA`
 *  * `GLFW_EXPOSE_NATIVE_X11`
 *  * `GLFW_EXPOSE_NATIVE_WAYLAND`
 *  * `GLFW_EXPOSE_NATIVE_MIR`
 *
 *  The available context API macros are:
 *  * `GLFW_EXPOSE_NATIVE_WGL`
 *  * `GLFW_EXPOSE_NATIVE_NSGL`
 *  * `GLFW_EXPOSE_NATIVE_GLX`
 *  * `GLFW_EXPOSE_NATIVE_EGL`
 *
 *  These macros select which of the native access functions that are declared
 *  and which platform-specific headers to include.  It is then up your (by
 *  definition platform-specific) code to handle which of these should be
 *  defined.
 */

/*************************************************************************
 * System headers and types
 *************************************************************************/

#if defined(GLFW_EXPOSE_NATIVE_WIN32)
 // This is a workaround for the fact that glfw3.h needs to export APIENTRY (for
 // example to allow applications to correctly declare a GL_ARB_debug_output
 // callback) but windows.h assumes no one will define APIENTRY before it does
 #undef APIENTRY
 #include <windows.h>
#elif defined(GLFW_EXPOSE_NATIVE_COCOA)
 #include <ApplicationServices/ApplicationServices.h>
 #if defined(__OBJC__)
  #import <Cocoa/Cocoa.h>
 #else
  typedef void* id;
 #endif
#elif defined(GLFW_EXPOSE_NATIVE_X11)
 #include <X11/Xlib.h>
 #include <X11/extensions/Xrandr.h>
#elif defined(GLFW_EXPOSE_NATIVE_WAYLAND)
 #include <wayland-client.h>
#elif defined(GLFW_EXPOSE_NATIVE_MIR)
 #include <mir_toolkit/mir_client_library.h>
#endif

#if defined(GLFW_EXPOSE_NATIVE_WGL)
 /* WGL is declared by windows.h */
#endif
#if defined(GLFW_EXPOSE_NATIVE_NSGL)
 /* NSGL is declared by Cocoa.h */
#endif
#if defined(GLFW_EXPOSE_NATIVE_GLX)
 #include <GL/glx.h>
#endif
#if defined(GLFW_EXPOSE_NATIVE_EGL)
 #include <EGL/egl.h>
#endif

/*************************************************************************
 * Functions
 *************************************************************************/

#if defined(GLFW_EXPOSE_NATIVE_WIN32)
/*! @brief Returns the adapter device name of the specified monitor.
 *
 *  @return The UTF-8 encoded adapter device name (for example `\\.\DISPLAY1`)
 *  of the specified monitor, or `NULL` if an [error](@ref error_handling)
 *  occurred.
 *
 *  @thread_safety This function may be called from any thread.  Access is not
 *  synchronized.
 *
 *  @since Added in version 3.1.
 *
 *  @ingroup native
 */
GLFWAPI const char* glfwGetWin32Adapter(GLFWmonitor* monitor);

/*! @brief Returns the display device name of the specified monitor.
 *
 *  @return The UTF-8 encoded display device name (for example
 *  `\\.\DISPLAY1\Monitor0`) of the specified monitor, or `NULL` if an
 *  [error](@ref error_handling) occurred.
 *
 *  @thread_safety This function may be called from any thread.  Access is not
 *  synchronized.
 *
 *  @since Added in version 3.1.
 *
 *  @ingroup native
 */
GLFWAPI const char* glfwGetWin32Monitor(GLFWmonitor* monitor);

/*! @brief Returns the `HWND` of the specified window.
 *
 *  @return The `HWND` of the specified window, or `NULL` if an
 *  [error](@ref error_handling) occurred.
 *
 *  @thread_safety This function may be called from any thread.  Access is not
 *  synchronized.
 *
 *  @since Added in version 3.0.
 *
 *  @ingroup native
 */
GLFWAPI HWND glfwGetWin32Window(GLFWwindow* window);
#endif

#if defined(GLFW_EXPOSE_NATIVE_WGL)
/*! @brief Returns the `HGLRC` of the specified window.
 *
 *  @return The `HGLRC` of the specified window, or `NULL` if an
 *  [error](@ref error_handling) occurred.
 *
 *  @thread_safety This function may be called from any thread.  Access is not
 *  synchronized.
 *
 *  @since Added in version 3.0.
 *
 *  @ingroup native
 */
GLFWAPI HGLRC glfwGetWGLContext(GLFWwindow* window);
#endif

#if defined(GLFW_EXPOSE_NATIVE_COCOA)
/*! @brief Returns the `CGDirectDisplayID` of the specified monitor.
 *
 *  @return The `CGDirectDisplayID` of the specified monitor, or
 *  `kCGNullDirectDisplay` if an [error](@ref error_handling) occurred.
 *
 *  @thread_safety This function may be called from any thread.  Access is not
 *  synchronized.
 *
 *  @since Added in version 3.1.
 *
 *  @ingroup native
 */
GLFWAPI CGDirectDisplayID glfwGetCocoaMonitor(GLFWmonitor* monitor);

/*! @brief Returns the `NSWindow` of the specified window.
 *
 *  @return The `NSWindow` of the specified window, or `nil` if an
 *  [error](@ref error_handling) occurred.
 *
 *  @thread_safety This function may be called from any thread.  Access is not
 *  synchronized.
 *
 *  @since Added in version 3.0.
 *
 *  @ingroup native
 */
GLFWAPI id glfwGetCocoaWindow(GLFWwindow* window);
#endif

#if defined(GLFW_EXPOSE_NATIVE_NSGL)
/*! @brief Returns the `NSOpenGLContext` of the specified window.
 *
 *  @return The `NSOpenGLContext` of the specified window, or `nil` if an
 *  [error](@ref error_handling) occurred.
 *
 *  @thread_safety This function may be called from any thread.  Access is not
 *  synchronized.
 *
 *  @since Added in version 3.0.
 *
 *  @ingroup native
 */
GLFWAPI id glfwGetNSGLContext(GLFWwindow* window);
#endif

#if defined(GLFW_EXPOSE_NATIVE_X11)
/*! @brief Returns the `Display` used by GLFW.
 *
 *  @return The `Display` used by GLFW, or `NULL` if an
 *  [error](@ref error_handling) occurred.
 *
 *  @thread_safety This function may be called from any thread.  Access is not
 *  synchronized.
 *
 *  @since Added in version 3.0.
 *
 *  @ingroup native
 */
GLFWAPI Display* glfwGetX11Display(void);

/*! @brief Returns the `RRCrtc` of the specified monitor.
 *
 *  @return The `RRCrtc` of the specified monitor, or `None` if an
 *  [error](@ref error_handling) occurred.
 *
 *  @thread_safety This function may be called from any thread.  Access is not
 *  synchronized.
 *
 *  @since Added in version 3.1.
 *
 *  @ingroup native
 */
GLFWAPI RRCrtc glfwGetX11Adapter(GLFWmonitor* monitor);

/*! @brief Returns the `RROutput` of the specified monitor.
 *
 *  @return The `RROutput` of the specified monitor, or `None` if an
 *  [error](@ref error_handling) occurred.
 *
 *  @thread_safety This function may be called from any thread.  Access is not
 *  synchronized.
 *
 *  @since Added in version 3.1.
 *
 *  @ingroup native
 */
GLFWAPI RROutput glfwGetX11Monitor(GLFWmonitor* monitor);

/*! @brief Returns the `Window` of the specified window.
 *
 *  @return The `Window` of the specified window, or `None` if an
 *  [error](@ref error_handling) occurred.
 *
 *  @thread_safety This function may be called from any thread.  Access is not
 *  synchronized.
 *
 *  @since Added in version 3.0.
 *
 *  @ingroup native
 */
GLFWAPI Window glfwGetX11Window(GLFWwindow* window);
#endif

#if defined(GLFW_EXPOSE_NATIVE_GLX)
/*! @brief Returns the `GLXContext` of the specified window.
 *
 *  @return The `GLXContext` of the specified window, or `NULL` if an
 *  [error](@ref error_handling) occurred.
 *
 *  @thread_safety This function may be called from any thread.  Access is not
 *  synchronized.
 *
 *  @since Added in version 3.0.
 *
 *  @ingroup native
 */
GLFWAPI GLXContext glfwGetGLXContext(GLFWwindow* window);

/*! @brief Returns the `GLXWindow` of the specified window.
 *
 *  @return The `GLXWindow` of the specified window, or `None` if an
 *  [error](@ref error_handling) occurred.
 *
 *  @thread_safety This function may be called from any thread.  Access is not
 *  synchronized.
 *
 *  @since Added in version 3.2.
 *
 *  @ingroup native
 */
GLFWAPI GLXWindow glfwGetGLXWindow(GLFWwindow* window);
#endif

#if defined(GLFW_EXPOSE_NATIVE_WAYLAND)
/*! @brief Returns the `struct wl_display*` used by GLFW.
 *
 *  @return The `struct wl_display*` used by GLFW, or `NULL` if an
 *  [error](@ref error_handling) occurred.
 *
 *  @thread_safety This function may be called from any thread.  Access is not
 *  synchronized.
 *
 *  @since Added in version 3.2.
 *
 *  @ingroup native
 */
GLFWAPI struct wl_display* glfwGetWaylandDisplay(void);

/*! @brief Returns the `struct wl_output*` of the specified monitor.
 *
 *  @return The `struct wl_output*` of the specified monitor, or `NULL` if an
 *  [error](@ref error_handling) occurred.
 *
 *  @thread_safety This function may be called from any thread.  Access is not
 *  synchronized.
 *
 *  @since Added in version 3.2.
 *
 *  @ingroup native
 */
GLFWAPI struct wl_output* glfwGetWaylandMonitor(GLFWmonitor* monitor);

/*! @brief Returns the main `struct wl_surface*` of the specified window.
 *
 *  @return The main `struct wl_surface*` of the specified window, or `NULL` if
 *  an [error](@ref error_handling) occurred.
 *
 *  @thread_safety This function may be called from any thread.  Access is not
 *  synchronized.
 *
 *  @since Added in version 3.2.
 *
 *  @ingroup native
 */
GLFWAPI struct wl_surface* glfwGetWaylandWindow(GLFWwindow* window);
#endif

#if defined(GLFW_EXPOSE_NATIVE_MIR)
/*! @brief Returns the `MirConnection*` used by GLFW.
 *
 *  @return The `MirConnection*` used by GLFW, or `NULL` if an
 *  [error](@ref error_handling) occurred.
 *
 *  @thread_safety This function may be called from any thread.  Access is not
 *  synchronized.
 *
 *  @since Added in version 3.2.
 *
 *  @ingroup native
 */
GLFWAPI MirConnection* glfwGetMirDisplay(void);

/*! @brief Returns the Mir output ID of the specified monitor.
 *
 *  @return The Mir output ID of the specified monitor, or zero if an
 *  [error](@ref error_handling) occurred.
 *
 *  @thread_safety This function may be called from any thread.  Access is not
 *  synchronized.
 *
 *  @since Added in version 3.2.
 *
 *  @ingroup native
 */
GLFWAPI int glfwGetMirMonitor(GLFWmonitor* monitor);

/*! @brief Returns the `MirSurface*` of the specified window.
 *
 *  @return The `MirSurface*` of the specified window, or `NULL` if an
 *  [error](@ref error_handling) occurred.
 *
 *  @thread_safety This function may be called from any thread.  Access is not
 *  synchronized.
 *
 *  @since Added in version 3.2.
 *
 *  @ingroup native
 */
GLFWAPI MirSurface* glfwGetMirWindow(GLFWwindow* window);
#endif

#if defined(GLFW_EXPOSE_NATIVE_EGL)
/*! @brief Returns the `EGLDisplay` used by GLFW.
 *
 *  @return The `EGLDisplay` used by GLFW, or `EGL_NO_DISPLAY` if an
 *  [error](@ref error_handling) occurred.
 *
 *  @thread_safety This function may be called from any thread.  Access is not
 *  synchronized.
 *
 *  @since Added in version 3.0.
 *
 *  @ingroup native
 */
GLFWAPI EGLDisplay glfwGetEGLDisplay(void);

/*! @brief Returns the `EGLContext` of the specified window.
 *
 *  @return The `EGLContext` of the specified window, or `EGL_NO_CONTEXT` if an
 *  [error](@ref error_handling) occurred.
 *
 *  @thread_safety This function may be called from any thread.  Access is not
 *  synchronized.
 *
 *  @since Added in version 3.0.
 *
 *  @ingroup native
 */
GLFWAPI EGLContext glfwGetEGLContext(GLFWwindow* window);

/*! @brief Returns the `EGLSurface` of the specified window.
 *
 *  @return The `EGLSurface` of the specified window, or `EGL_NO_SURFACE` if an
 *  [error](@ref error_handling) occurred.
 *
 *  @thread_safety This function may be called from any thread.  Access is not
 *  synchronized.
 *
 *  @since Added in version 3.0.
 *
 *  @ingroup native
 */
GLFWAPI EGLSurface glfwGetEGLSurface(GLFWwindow* window);
#endif

#ifdef __cplusplus
}
#endif

#endif /* _glfw3_native_h_ */


text:
# Contributing

Your contributions are always welcome!

## Guidelines

- Add commits to Pull Request.

- Add the link: The format for the item is [software link] - [Paper link] (Available in which package) and the headings may contain a ink to a review paper.

- Add a section if needed.

  - Add the section title to Table of Contents.

- Check your spelling and grammar.

- Remove any trailing whitespace.

repo_name: phonchi/Computational-CryoEM

repo_language: Jupyter Notebook

mime_type: text/plain


text:
import numpy as np
import matplotlib.pyplot as plt
import numpy.linalg as LA
from scipy.sparse.linalg import svds, eigs
from numpy.linalg import svd
from sklearn import metrics
from datetime import datetime as dt
import os
import pandas as pd
import sys
import seaborn as sns

import mrc
from mrc import LazyImage
import mrcfile

### Some of code are from https://github.com/PrincetonUniversity/ASPIRE-Python and https://github.com/zhonge/cryodrgn

########### io #################
class Starfile():

    def __init__(self, headers, df):
        self.headers = headers
        self.df = df

    @classmethod
    def load(self, starfile):
        f = open(starfile,'r')
        # get to data block
        while 1:
            for line in f:
                if line.startswith('data_'):
                    break
            break
        # get to header loop
        while 1:
            for line in f:
                if line.startswith('loop_'):
                    break
            break
        # get list of column headers
        while 1:
            headers = []
            for line in f:
                if line.startswith('_'):
                    headers.append(line)
                else:
                    break
            break
        # assume the rest is the body
        headers = [h.strip().split()[0] for h in headers]
        body = [line]
        for lines in f:
            body.append(lines)
        # remove last line of body if empty
        if body[-1].strip() == '':
            body = body[:-1]
        # put data into an array and instantiate as dataframe
        words = [l.strip().split() for l in body]
        words = np.array(words)
        data = {h:words[:,i] for i,h in enumerate(headers)}
        df = pd.DataFrame(data=data)
        return self(headers, df)

    def write(self, outstar):
        f = open(outstar,'w')
        f.write('# Created {}\n'.format(dt.now()))
        f.write('\n')
        f.write('data_\n\n')
        f.write('loop_\n')
        self.headers = list(self.df)
        for idx, h in enumerate(self.headers):
            f.write(h+' #%s'%(idx+1))
            f.write('\n')
        for i in self.df.index:
            f.write(' '.join(self.df.loc[i]))
            f.write('\n')
        #f.write('\n'.join([' '.join(self.df.loc[i]) for i in range(len(self.df))]))

    def get_particles(self, datadir=None, lazy=False):
        '''
        Return particles of the starfile

        Input:
            datadir (str): Overwrite base directories of particle .mrcs
                Tries both substituting the base path and prepending to the path
            If lazy=True, returns list of LazyImage instances, else np.array
        '''
        particles = self.df['_rlnImageName']

        # format is index@path_to_mrc
        particles = [x.split('@') for x in particles]
        ind = [int(x[0])-1 for x in particles] # convert to 0-based indexing
        mrcs = [x[1] for x in particles]
        if datadir is not None:
            mrcs = prefix_paths(mrcs, datadir)
        for path in set(mrcs):
            assert os.path.exists(path), "%s not found"%path
        D = mrc.parse_header(mrcs[0]).D # image size along one dimension in pixels
        dtype = np.float32
        stride = np.float32().itemsize*D*D
        dataset = [LazyImage(f, (D,D), dtype, 1024+ii*stride) for ii,f in zip(ind, mrcs)]
        if not lazy:
            dataset = np.array([x.get() for x in dataset])
        return dataset

def prefix_paths(mrcs, datadir):
    mrcs1 = ['{}/{}'.format(datadir, os.path.basename(x)) for x in mrcs]
    mrcs2 = ['{}/{}'.format(datadir, x) for x in mrcs]
    try:
        for path in set(mrcs1):
            assert os.path.exists(path)
        mrcs = mrcs1
    except:
        for path in set(mrcs2):
            assert os.path.exists(path), "%s not found"%path
        mrcs = mrcs2
    return mrcs

def csparc_get_particles(csfile, datadir=None, lazy=True):
    metadata = np.load(csfile)
    ind = metadata['blob/idx'] # 0-based indexing
    mrcs = metadata['blob/path'].astype(str).tolist()
    if datadir is not None:
        mrcs = prefix_paths(mrcs, datadir)
    for path in set(mrcs):
        assert os.path.exists(path), "%s not found"%path
    D = metadata[0]['blob/shape'][0]
    dtype = np.float32
    stride = np.float32().itemsize*D*D
    dataset = [LazyImage(f, (D,D), dtype, 1024+ii*stride) for ii,f in zip(ind, mrcs)]
    if not lazy:
        dataset = np.array([x.get() for x in dataset])
    return dataset

########### util #################
def log(msg):
    print('{}     {}'.format(dt.now().strftime('%Y-%m-%d %H:%M:%S'), msg))
    sys.stdout.flush()

def print_ctf_params(params):
    assert len(params) == 9
    log('Image size (pix)  : {}'.format(int(params[0])))
    log('A/pix             : {}'.format(params[1]))
    log('DefocusU (A)      : {}'.format(params[2]))
    log('DefocusV (A)      : {}'.format(params[3]))
    log('Dfang (deg)       : {}'.format(params[4]))
    log('voltage (kV)      : {}'.format(params[5]))
    log('cs (mm)           : {}'.format(params[6]))
    log('w                 : {}'.format(params[7]))
    log('Phase shift (deg) : {}'.format(params[8]))

def parse_ctf_star(df, D, angpix=None):
    N = len(df)
    if angpix == None:
        if set(['_rlnDetectorPixelSize','_rlnMagnification']).issubset(df.columns):
            Apix = float(df['_rlnDetectorPixelSize'][0])*10000/float(df['_rlnMagnification'][0]);
        else:
            Apix = 1
    else:
        Apix = angpix
    ctf_params = np.zeros((N, 9))
    ctf_params[:,0] = D
    ctf_params[:,1] = Apix

    for i,header in enumerate(['_rlnDefocusU', '_rlnDefocusV', '_rlnDefocusAngle', '_rlnVoltage', '_rlnSphericalAberration', '_rlnAmplitudeContrast', '_rlnPhaseShift']):
        ctf_params[:,i+2] = df[header]
    log('CTF parameters for first particle:')
    print_ctf_params(ctf_params[0])
    return ctf_params

def parse_pose_star(df):
    # parse rotations
    N = len(df)
    keys = ('_rlnAngleRot','_rlnAngleTilt','_rlnAnglePsi')
    euler = np.empty((N,3))
    euler[:,0] = df['_rlnAngleRot']
    euler[:,1] = df['_rlnAngleTilt']
    euler[:,2] = df['_rlnAnglePsi']
    log('Euler angles (Rot, Tilt, Psi):')
    log(euler[0])
    log('Converting to rotation matrix:')
    rot = np.asarray([R_from_relion(*x) for x in euler])
    log(rot[0])

    # parse translations
    trans = np.empty((N,2))
    trans[:,0] = df['_rlnOriginX']
    trans[:,1] = df['_rlnOriginY']
    log('Translations:')
    log(trans[0])
    return (euler, trans, rot)

def R_from_relion(a,b,y):
    a *= np.pi/180.
    b *= np.pi/180.
    y *= np.pi/180.
    ca, sa = np.cos(a), np.sin(a)
    cb, sb = np.cos(b), np.sin(b)
    cy, sy = np.cos(y), np.sin(y)
    Ra = np.array([[ca,-sa,0],[sa,ca,0],[0,0,1]])
    Rb = np.array([[cb,0,-sb],[0,1,0],[sb,0,cb]])
    Ry = np.array(([cy,-sy,0],[sy,cy,0],[0,0,1]))
    R = np.dot(np.dot(Ry,Rb),Ra)
    R[0,1] *= -1
    R[1,0] *= -1
    R[1,2] *= -1
    R[2,1] *= -1
    return R

def visualise_images(X, n_images, n_columns, randomise=True):
    indices = np.arange(X.shape[0])
    if randomise:
        np.random.shuffle(indices)
    indices = indices[:n_images]
    cmap = plt.cm.Greys_r
    n_rows = np.ceil(n_images / n_columns)
    fig = plt.figure(figsize=(2*n_columns, 2*n_rows))
    fig.subplots_adjust(left=0, right=1, bottom=0, top=1, hspace=0.05, wspace=0.05)

    # plot the digits: each image is 8x8 pixels
    for i, e in enumerate(indices):
        ax = fig.add_subplot(n_rows, n_columns, i + 1, xticks=[], yticks=[])
        ax.imshow(X[e], cmap=cmap, interpolation='nearest')

def matlab2py(i_matrix):
    tmp = np.swapaxes(i_matrix,0,2)
    return np.swapaxes(tmp,1,2).copy()

def save_mrc(images, filename):
    """Shortcut for writing out MRC files."""
    o = mrcfile.new(filename, overwrite=True)
    o.set_data(images.astype(np.float32))
    o.close()

def cart2pol(x, y):
    """
    Convert Cartesian to Polar Coordinates. All input arguments must be the same shape.
    :param x: x-coordinate in Cartesian space
    :param y: y-coordinate in Cartesian space
    :return: A 2-tuple of values:
        theta: angular coordinate/azimuth
        r: radial distance from origin
    """
    return np.arctan2(y, x), np.hypot(x, y)

def grid_2d(n, shifted=False, normalized=True):
    """
    Generate two dimensional grid.

    :param n: the number of grid points in each dimension.
    :param shifted: shifted by half of grid or not when n is even.
    :param normalized: normalize the grid in the range of (-1, 1) or not.
    :return: the rectangular and polar coordinates of all grid points.
    """
    grid = np.ceil(np.arange(-n/2, n/2))

    if shifted and n % 2 == 0:
        grid = np.arange(-n/2+1/2, n/2+1/2)

    if normalized:
        if shifted and n % 2 == 0:
            grid = grid / (n/2-1/2)
        else:
            grid = grid / (n/2)

    x, y = np.meshgrid(grid, grid, indexing='ij')
    phi, r = cart2pol(x, y)

    return {
        'x': x,
        'y': y,
        'phi': phi,
        'r': r
    }

def estimate_noise(in_matrix):
    ### in_matrix in matlab format
    arr4 = np.swapaxes(in_matrix,0,2)
    arr4 = np.swapaxes(arr4,0,1)
    g2d = grid_2d(in_matrix.shape[1])
    mask = g2d['r'] >= 1

    first_moment = 0
    second_moment = 0

    images_masked = arr4 * np.expand_dims(mask, 2)

    _denominator = in_matrix.shape[0] * np.sum(mask)
    first_moment += np.sum(images_masked) / _denominator
    second_moment += np.sum(np.abs(images_masked**2)) / _denominator
    return second_moment - first_moment**2

############ analysis #################
def plot_euler(euler,trans,plot_psi=True,plot_trans=True):
    theta = euler[:,0]
    phi = euler[:,1]
    psi = euler[:,2]
    hexplot = sns.jointplot(theta,phi,kind='hex',
              xlim=(-180,180),
              ylim=(0,180)).set_axis_labels("theta", "phi")
    plt.subplots_adjust(left=0.1, right=0.9, top=0.9, bottom=0.1)
    cbar = hexplot.fig.add_axes([.9,.1,.04, .7])
    plt.colorbar(cax=cbar)
    plt.show()
    if plot_psi:
        plt.figure()
        plt.hist(psi)
        plt.xlabel('psi')
    if plot_trans:
        sns.jointplot(trans[:,0],trans[:,1],
              kind='hex').set_axis_labels('tx','ty')

def plot_defocus(ctfs):
    plt.hist(ctfs[:,2])
    plt.xlabel('DefocusU (um)')
    plt.figure()
    plt.hist(ctfs[:,3])
    plt.xlabel('DefocusV (um)')

def compute_ctf_np(freqs, dfu, dfv, dfang, volt, cs, w, phase_shift=0, bfactor=None):
    '''
    Compute the 2D CTF

    Input:
        freqs (np.ndarray) Nx2 array of 2D spatial frequencies
        dfu (float): DefocusU (Angstrom)
        dfv (float): DefocusV (Angstrom)
        dfang (float): DefocusAngle (degrees)
        volt (float): accelerating voltage (kV)
        cs (float): spherical aberration (mm)
        w (float): amplitude contrast ratio
        phase_shift (float): degrees
        bfactor (float): envelope fcn B-factor (Angstrom^2)
    '''
    # convert units
    volt = volt * 1000
    cs = cs * 10**7
    dfang = dfang * np.pi / 180
    phase_shift = phase_shift * np.pi / 180

    # lam = sqrt(h^2/(2*m*e*Vr)); Vr = V + (e/(2*m*c^2))*V^2
    lam = 12.2639 / np.sqrt(volt + 0.97845e-6 * volt**2)
    x = freqs[:,0]
    y = freqs[:,1]
    ang = np.arctan2(y,x)
    s2 = x**2 + y**2
    df = .5*(dfu + dfv + (dfu-dfv)*np.cos(2*(ang-dfang)))
    gamma = 2*np.pi*(-.5*df*lam*s2 + .25*cs*lam**3*s2**2) - phase_shift
    ctf = np.sqrt(1-w**2)*np.sin(gamma) - w*np.cos(gamma)
    if bfactor is not None:
        ctf *= np.exp(-bfactor/4*s2)
    return np.require(ctf,dtype=freqs.dtype)

def plot_ctf(ctf_params):
    assert len(ctf_params) == 9
    import matplotlib.pyplot as plt
    import seaborn as sns
    D = int(ctf_params[0])
    Apix = ctf_params[1]
    freqs = np.stack(np.meshgrid(np.linspace(-.5,.5,D,endpoint=False),np.linspace(-.5,.5,D,endpoint=False)),-1)/Apix
    freqs = freqs.reshape(-1,2)
    c = compute_ctf_np(freqs, *ctf_params[2:])
    sns.heatmap(c.reshape(D, D))

############ metrics#################

def purity_score(y_true, y_pred):
    # compute contingency matrix (also called confusion matrix)
    contingency_matrix = metrics.cluster.contingency_matrix(y_true, y_pred)
    # return purity
    return np.sum(np.amax(contingency_matrix, axis=0)) / np.sum(contingency_matrix)

def c_purity_score(y_true, y_pred):
    # compute contingency matrix (also called confusion matrix)
    contingency_matrix = metrics.cluster.contingency_matrix(y_true, y_pred)
    # return purity
    return np.sum(np.amax(contingency_matrix, axis=1)) / np.sum(contingency_matrix)


text:
# Computational CryoEM Methods
[![Awesome](https://awesome.re/badge.svg)](https://awesome.re)

This repository is a curated list of computational cryo-EM methods mainly targeted for single particle analysis! If you are looking for computational methods for cryo-ET, see [here](https://github.com/phonchi/Computational-CryoET).

You will find the paper and associated software for the popular algorithms used in the cryo-EM field.

The format for the item is [software link] - [Paper link] (Available in which package) and the headings may contain a link to a review paper.

- [**Computational CryoEM Methods**](#computational-cryoEM-methods)
- [**Resources**](#resources)
- [**Workflow**](#workflow)
    - [Motion correction](#motion-correction)
    - [CTF estimation](#ctf-estimation)
    - [Particle picking](#particle-picking)
    - [2D classification](#2D-classification)
    - [Ab-initial Model](#ab-initial-model)
    - [3D refinement](#3d-refinement)
    - [3D variability analysis](#3D-variability-analysis)
    - [Postprocessing](#postprocessing)
- [**Conventions**](#conventions)
- [**Tips**](#tips)
    - [Format Conversion](#format-conversion)
    - [Performing Focus Classification](#performing-focus-classification)
    - [Display images](#display-images)
    - [Using Chimera](#using-chimera)
- [**Contributing**](#contributing)

---

# Resources
## Preface
* [Protein Folding Problem](https://www.youtube.com/watch?v=KpedmJdrTpY&t=2s&ab_channel=DeepMind) - Great video describes the importance of understanding protein structure.
* [What is cryoEM](https://www.youtube.com/watch?v=Qq8DO-4BnIY&ab_channel=UCSanFrancisco%28UCSF%29) - Great video that describes how cryoEM can help us understand the protein structure.

## Introduction
* [A quick introduction video for cryo-EM](https://www.youtube.com/watch?v=Qq8DO-4BnIY)
* [Cryo-EM 101](https://cryoem101.org/) - Deomnstrate the principles of cryo-EM using a media-rich approach with videos, animations, interactive simulation.
* [Cryo-EM workflow by Nvidia](https://resources.nvidia.com/en-us-drug-discovery/) - Describe the blueprint of cryo-EM using GPU processing.
* [A Primer to Single-Particle Cryo-Electron Microscopy](https://www.sciencedirect.com/science/article/pii/S0092867415003700?via%3Dihub) - A great review to start with.
* [Principles of cryo-EM single-particle image processing](https://pubmed.ncbi.nlm.nih.gov/26705325/) - Highlight the principles when dealing with a cryo-EM dataset.
* [Single-particle cryo-electron microscopy: Mathematical theory, computational challenges, and opportunities](https://web.math.princeton.edu/~amits/publications/IEEE_Signal_Processing_Magazine.pdf) - Review that covers recent updates.
* [Computational Methods for Single-Particle Electron Cryomicroscopy](https://www.annualreviews.org/doi/abs/10.1146/annurev-biodatasci-021020-093826) - A comprehensive review for the mathematical background behind single-particle cryo-EM methods.
* [Cryo-EM: Breakthroughs in Chemistry, Structural Biology, and Statistical Underpinnings](https://statistics.stanford.edu/research/cryo-em-breakthroughs-chemistry-structural-biology-and-statistical-underpinnings) - Review that covers statistical problems behind cryo-EM.

## Courses
* [Cryo-EM Principles](https://cryoemprinciples.yale.edu/)
* [Getting Started in Cryo-EM](http://cryo-em-course.caltech.edu/overview)
* [Electron Microscopy Training by MRC LAB](https://www2.mrc-lmb.cam.ac.uk/research/scientific-training/electron-microscopy/)
* [nysbc training course](http://semc.nysbc.org/the-winter-spring-2018-em-course/)

## Resources
* [3DEM Methods](http://3demmethods.i2pc.es/index.php/Main_Page) - A great wiki that collects papers or books for computational methods.
* [A collective resource](https://github.com/barrykui/awesome-cryoem) - A great repository that covers single particle analysis, model building and tomography.
* [Math behind CryoEM](https://github.com/geoffwoollard/learn_cryoem_math) - A great repository that collects the materials which elaborate the math behind single-particle analysis.
* [Chimera Tutorial](http://embo2015.cryst.bbk.ac.uk/embo2015/course/practicals/prac-7/Practical-7-Maya-Topf-9-Sept-15.pdf)

## Events and News
* [CryoEM and AI reveal a structure of SARS-CoV-2 Nsp2, a multifunctional protein involved in key host processes](https://www.biorxiv.org/content/10.1101/2021.05.10.443524v1)
* [3DEM Events](http://www.emdatabank.org/3dem_events.html)
* [NRAMM Events](http://nramm.nysbc.org/workshops-and-courses/?lcp_page0=1#lcp_instance_0)

## [Software](http://3dem.ucsd.edu/software.shtm)
* [Scipion](http://scipion.i2pc.es/) - An integrated platform that allows users to use a variety of methods in the same framework. [[Documentation](https://scipion-em.github.io/docs/)]
* [Relion](https://www3.mrc-lmb.cam.ac.uk/relion/index.php?title=Main_Page) - A comprehensive package that utilizes Bayesian approach for 2D classification and 3D refinement. [[Documentation](https://relion.readthedocs.io/en/latest/)]
* [CryoSparc 2](https://cryosparc.com/) - A general package that employs stochastic gradient descent, branch and bound as well as GPU acceleration for rapid reconstructions.
* [Sphire](http://sphire.mpg.de/) - A general package that contains neural network methods for particle picking, denoising and classification.
* [CisTEM](https://github.com/timothygrant80/cisTEM) - An easy-to-use framework that implements a complete pipeline for single-particle analysis.
* [EMAN2](https://blake.bcm.edu/emanwiki/EMAN2) - A comprehensive package that contains a python interface and handy scripts for common tasks.
* [SPIDER](https://spider.wadsworth.org/spider_doc/spider/docs/spider.html) - A well-known package that implements image processing methods for electron microscopy.
* [ASPIRE](http://spr.math.princeton.edu/) - A software that uses algorithms based on rigorous mathematical theory.

# Workflow
## Import Data
1. From [EMPIAR](https://www.ebi.ac.uk/pdbe/emdb/empiar/) or [EMDB](https://www.ebi.ac.uk/pdbe/emdb/index.html/).

 DataSet (Molecule) | File Size | Movie (#frame) | Micrograph Size (Pixel Size) | Picked Particles (Size) | Final Resolution(A)(Symmetry)| Notes
 :--------: | :--------: | :--------: | :--------: | :--------: | :--------: | :--------:
 [70S ribosome](https://www3.mrc-lmb.cam.ac.uk/relion/index.php?title=Classification_example)  | 0.7GB     | N/A     |  N/A (2.82)    |  10,000 (130x130)    |  ~9 (C1) | Test set for classification (2 classes)
 [50S ribosome](https://www.ebi.ac.uk/pdbe/emdb/empiar/entry/10076/)   | 50.3GB     | N/A     |  N/A (1.31)    |  131,899 (320X320)    | N/A| Test set for classification (4~7 classes)
 [80S ribosome](https://www3.mrc-lmb.cam.ac.uk/relion/index.php?title=Benchmarks_%26_computer_hardware)   | 1.2TB     | 1081(16)     |  4096x4096 (1.34)    | 105,247 (360x360)     | 3.2 (C1) | Test set for computational performance
 [Beta-Gal](https://www.ebi.ac.uk/pdbe/emdb/empiar/entry/10204/)   | 321.4GB   | 1338(49)     |  3710x3838 (0.885)    | N/A     | 2.6 (D2)| Test set for high-resolution reconstruction
 [Apoferritin](https://www.ebi.ac.uk/pdbe/emdb/empiar/entry/10200/)   | 191.5GB     | 1255(40)     |  3710x3838 (0.814)    |   N/A   | 1.65 (O)| Test set for high-resolution reconstruction
 [T20S](https://www.ebi.ac.uk/pdbe/emdb/empiar/entry/10025/)   | 2.0TB     | 196(38)     |  7420x7676 (0.6575)    |  49,955(448X448)    | 2.8 (D7)| Test set for high-resolution reconstruction
 [TRPV1](https://www.ebi.ac.uk/pdbe/emdb/empiar/entry/10059/)   | 93.8GB     | 1200(1)     |  3710x3838 (1.256)    |  218,805 (192X192)    | 2.95 (C4)| Test set for Membrane protein (Nonuniform reconstruction)
 [Spliceosome](https://www.ebi.ac.uk/pdbe/emdb/empiar/entry/10180/)   | 126.5GB     | N/A     |  N/A (1.699)    |  327,490 (320X320)    | N/A | Test set for continuous conformation

3. From tutorial data set.

 DataSet (Molecule) | File Size | Movie (#frame) | Micrograph Size (Pixel Size) | Picked Particles (Size) | Final Resolution(A)(Symmetry)
 :--------: | :--------: | :--------: | :--------: | :--------: | :--------:
 [Scipion](https://github.com/I2PC/scipion/wiki/tutorials/scipion_tutorial_betagal.pdf)
  ([Beta-gal](https://scipion.cnb.csic.es/downloads/scipion/data/tests/relion13_tutorial/))    | 4.0GB     | 15(16)     |  1950x1950 (3.54)    |  (100x100)    | 7.3 (D2)
 [Relion 3](https://relion.readthedocs.io/en/latest/SPA_tutorial/Introduction.html)
  ([Beta-gal](https://www.ebi.ac.uk/pdbe/emdb/empiar/entry/10204/))    | 3.1GB     | 24(24)     |  3710x3838 (0.885)    | (256x256)     |2.9 (D2)
 [CryoSparc 2](https://cryosparc.com/docs/tutorials/t20s) ([T20S](https://www.ebi.ac.uk/pdbe/emdb/empiar/entry/10025/))    | 8.1GB     | 20(38)     |  7420x7676 (0.6575)    | (440x440)     | 2.93 (D7)
 [CisTEM](https://cistem.org/documentation#tab-1-1) ([ApoFerritin](https://www.ebi.ac.uk/pdbe/emdb/empiar/entry/10146/))    | 5.6GB     | 20(50)     |  1240x1200 (1.5)    |      | 3.0 (O)

3. Generate simulation data
    * [Discrete](Simulation_data.ipynb)
    * [Continuous](https://github.com/evanseitz/cryoEM_synthetic_continua)
    * From [CryoDRGN](https://github.com/zhonge/cryodrgn) - Use [`write_star.py`](https://github.com/zhonge/cryodrgn/blob/master/cryodrgn/commands_utils/write_star.py) to convert to star format

DataSet (Molecule) | File Size | Micrograph Size (Pixel Size) | Picked Particles (Size)
 :--------: | :--------: | :--------: | :--------:
 [Ribosome](https://zenodo.org/record/4355284#.YHPD3a8zbmF)  | 3.1GB     | N/A (3)   | 50,000 (128x128)
 [Uniform](https://zenodo.org/record/4355284#.YHPD3a8zbmF)   | 3.0GB     | N/A (6)   | 50,000 (128x128)
 [Cooperative](https://zenodo.org/record/4355284#.YHPD3a8zbmF)  | 3.0GB     | N/A (6)   | 50,000 (128x128)
 [Noncontiguous](https://zenodo.org/record/4355284#.YHPD3a8zbmF)   | 3.0GB     | N/A (6)   | 50,000 (128x128)

## Image Formation model
1. [`TEM Simulator`](http://tem-simulator.sourceforge.net/) - [**Simulation of transmission electron microscope images of biological specimens**](https://onlinelibrary.wiley.com/doi/10.1111/j.1365-2818.2011.03497.x)

2. [`InSilicoTEM`](https://github.com/M4I-nanoscopy/InSilicoTEM/tree/v2.0.0) - [**Image formation modeling in cryo-electron microscopy**](https://www.sciencedirect.com/science/article/abs/pii/S1047847713001226)

3. [`CisTEM_Simulate`](https://github.com/timothygrant80/cisTEM/blob/master/src/programs/simulate/simulate.cpp) - [**Cryo-TEM simulations of amorphous radiation-sensitive samples using multislice wave propagation**](https://www.biorxiv.org/content/10.1101/2021.02.19.431636v2)

## [Motion correction](https://www.sciencedirect.com/science/article/pii/S0076687916300271?via%3Dihub)

* ### Whole frame alignment
1. [`Unblur`](https://grigoriefflab.umassmed.edu/unblur_summovie) - [**Measuring the optimal exposure for single particle cryo-EM using a 2.6 Å reconstruction of rotavirus VP6**](https://elifesciences.org/articles/06980) ( Use in [`CisTEM`](https://cistem.org/documentation#tab-1-6))
2. `Full-frame motion correction` - (Use in `CryoSparc`)
3. `Optical Flow` - [**Alignment of direct detection device micrographs using a robust Optical Flow approach**](https://www.sciencedirect.com/science/article/pii/S1047847715000313) (Use in [`Xmipp`](http://xmipp.i2pc.es/))

* ### Patch-based/Per-particle based alignment
1. [`MotionCorr2`](https://emcore.ucsf.edu/ucsf-motioncor2) - [**MotionCor2: anisotropic correction of beam-induced motion for improved cryo-electron microscopy**](https://www.nature.com/articles/nmeth.4193?proof=t) (Use in `Relion` and `CryoSparc`)
2. [`Alignparts`](https://sites.google.com/site/rubinsteingroup/direct-detector-align_lmbfgs) - [**Alignment of cryo-EM Movies of Individual Particles by Optimization of Image Translations**](https://pubmed.ncbi.nlm.nih.gov/26296328/) (Use in [`CryoSparc`(local, patch)](https://cryosparc.com/docs/tutorials/patch-motion-ctf))
3. [`Warp`](http://www.warpem.com/warp/?page_id=185) - Contains [Patch-based motion correction](http://www.warpem.com/warp/).
4. `FlexAlign` - [**FlexAlign: An Accurate and Fast Algorithm for Movie Alignment in Cryo-Electron Microscopy**](https://www.mdpi.com/2079-9292/9/6/1040/htm) - (Use in `Scipion`)

* ### [Damage compensation](https://www.sciencedirect.com/science/article/pii/S0076687910810158?via%3Dihub)
    * This functionality is available in most of the motion correction tool
    * The [**dark/gain corrected**](https://pubmed.ncbi.nlm.nih.gov/27572724/) is also conducted in this stage

## [CTF estimation](https://www.sciencedirect.com/science/article/pii/S0076687910820026?via%3Dihub)
* ### Whole frame
1. [`CTFFIND4`](https://grigoriefflab.umassmed.edu/ctffind4) - [**CTFFIND4: Fast and accurate defocus estimation from electron micrographs**](https://www.sciencedirect.com/science/article/pii/S1047847715300460) ( Use in [`CisTEM`](http://grigoriefflab.janelia.org/ctffind4))
2. [`gCTF`](https://www2.mrc-lmb.cam.ac.uk/research/locally-developed-software/zhang-software/#gctf) - [**Gctf: Real-time CTF determination and correction**](https://www.ncbi.nlm.nih.gov/pmc/articles/PMC4711343/)
* ### Patch-based/tilt data
1. `Patch-Based CTF Estimation` - [**Real-time cryo-electron microscopy data preprocessing with Warp**](https://www.nature.com/articles/s41592-019-0580-y) ( Use in [`CryoSparc`](https://cryosparc.com/docs/tutorials/patch-motion-ctf), `Warp`)
2. [`goCTF`](https://www.lsi.umich.edu/science/centers-technologies/cryo-electron-microscopy/research/goctf) - [**goCTF: Geometrically optimized CTF determination for single-particle cryo-EM**](https://www.lsi.umich.edu/science/centers-technologies/cryo-electron-microscopy/research/goctf)
3. [`Warp`](http://www.warpem.com/warp/?page_id=185) - Contains Patch-based CTF estimation.
4. [`novaCTF`](https://github.com/turonova/novaCTF) - [**Efficient 3D-CTF correction for cryo-electron tomography using NovaCTF improves subtomogram averaging resolution to 3.4 Å**](https://www.sciencedirect.com/science/article/pii/S1047847717301272?via%3Dihub)

* ### Denoising micrograph
1. [`Restore`](https://github.com/eugenepalovcak/restore) - [**Enhancing SNR and generating contrast for cryo-EM images with convolutional neural networks**](https://journals.iucr.org/m/issues/2020/06/00/pw5015/)
2. `Topaz` - [**Topaz-Denoise: general deep denoising models for cryoEM**](https://www.nature.com/articles/s41467-020-18952-1)
3.  [`JANNI`](https://sphire.mpg.de/wiki/doku.php?id=janni) - Just Another Noise 2 Noise Implementation
4.  `Warp` - Contains methods base on noise2noise and deconvolution filter

## Particle picking
* ### Semi-supervised picking
1. [`Topaz`](https://github.com/tbepler/topaz) - [**Positive-unlabeled convolutional neural networks for particle picking in cryo-electron micrographs**](https://www.nature.com/articles/s41592-019-0575-8). [[Video]](https://www.youtube.com/watch?v=FnCnsT3GIC4&list=PLFEB3YHuxu11Jp_pOCIEtXxSqozFHve0O&index=10)
2. [`Cryolo`](https://cryolo.readthedocs.io/en/stable/) - [**SPHIRE-crYOLO is a fast and accurate fully automated particle picker for cryo-EM**](https://www.nature.com/articles/s42003-019-0437-z). [[Video](https://www.youtube.com/watch?v=JTgldM4wAAk&list=UUbo1TnKiXGtkE_R5b526JmQ&index=6&t=1s&app=desktop)]
4. `Xmipp` - [**A pattern matching approach to the automatic selection of particles from low-contrast electron micrographs**](https://academic.oup.com/bioinformatics/article/29/19/2460/189951)
5. [`EPicker`](https://github.com/thuem/EPicker) - [**EPicker is an exemplar-based continual learning approach for knowledge accumulation in cryoEM particle picking**](https://www.nature.com/articles/s41467-022-29994-y)

* ### Template-based picking
1. [`Relion`](https://www.sciencedirect.com/science/article/pii/S1047847714002615), `CryoSparc` - Use 2D class averages or 3D projection for more accurate particle picking
* ### Automatic picking
1. [`DoG`](https://emg.nysbc.org/redmine/projects/software/wiki/DoGpicker) - [**DoG Picker and TiltPicker: software tools to facilitate particle selection in single particle electron microscopy**](https://www.ncbi.nlm.nih.gov/pmc/articles/PMC2768396/)
2. `LoG` - [**The Laplacian of Gaussian and Arbitrary Z-Crossings Approach Applied to Automated Single Particle Reconstruction**](https://pubmed.ncbi.nlm.nih.gov/17490891/) (Use in `Relion` auto-picking)
3. [`APPLE`](https://github.com/PrincetonUniversity/APPLEpicker) - [**APPLE picker: Automatic particle picking, a low-effort cryo-EM framework**](https://www.sciencedirect.com/science/article/pii/S1047847718302326) (Use in `ASPIRE` auto-picking)
4. [`KLT picker`](https://github.com/ShkolniskyLab/kltpicker) - [**KLT picker: Particle picking using data-driven optimal templates**](https://www.sciencedirect.com/science/article/pii/S1047847720300356)

* ### Denoising particle and dimenison reduction
1. [`2SDR`](http://sabid.stat.sinica.edu.tw/start) - [**Two-stage dimension reduction for noisy high-dimensional images and application to Cryogenic Electron Microscopy**](https://www.intlpress.com/site/pub/pages/journals/items/amsa/content/vols/0005/0002/a004/index.php)
2. [`CWF`](https://github.com/PrincetonUniversity/cwf_denoise) - [**Denoising and covariance estimation of single particle cryo-EM images**](https://www.sciencedirect.com/science/article/pii/S104784771630082X) (Use in `ASPIRE`)
3. [`GAN`](https://github.com/cianfrocco-lab/GAN-for-Cryo-EM-image-denoising) - [**Generative adversarial networks as a tool to recover structural information from cryo-electron microscopy data**](https://www.biorxiv.org/content/10.1101/256792v1.article-info)

## 2D classification

* ### Multirefence alignment based classification
1. `ISAC` - [**Iterative Stable Alignment and Clustering of 2D Transmission Electron Microscope Images**](https://www.ncbi.nlm.nih.gov/pmc/articles/PMC3426367/). [[GPU version]](https://sphire.mpg.de/wiki/doku.php?id=gpu_isac)
2. `CL2D` - [**A clustering approach to multireference alignment of single-particle projections in electron microscopy**](https://www.ncbi.nlm.nih.gov/pmc/articles/PMC2893300/)

* ### [Maximum-likelihood based classification](https://www.sciencedirect.com/science/article/pii/S0076687910820117?via%3Dihub)
1. `Relion` - Bayesian (Empirical Bayes) approach
2. `CryoSparc` - Bayesian with Branch and bound method
3. `CisTEM` - Maximum likelihood method
4. [`SubspaceEM`](https://www.mathworks.com/matlabcentral/fileexchange/50091-subspaceem-a-fast-maximum-a-posteriori-algorithm-for-cryo-em-single-particle-reconstruction) - [**SubspaceEM: A fast maximum-a-posteriori algorithm for cryo-EM single particle reconstruction**](https://www.sciencedirect.com/science/article/pii/S1047847715000714)

* ### Mixed approach classification
6. [`ROME`](http://ipccsb.dfci.harvard.edu/rome/) - [**Massively parallel unsupervised single-particle cryo-EM data clustering via statistical manifold learning**](https://journals.plos.org/plosone/article?id=10.1371/journal.pone.0182130)

* ### Automatic selection of 2D classes
1. [`Cryoassess`](https://github.com/cianfrocco-lab/Automatic-cryoEM-preprocessing) -[**High-Throughput Cryo-EM Enabled by User-Free Preprocessing Routines**](https://www.sciencedirect.com/science/article/pii/S0969212620300800?via%3Dihub)
2. [`Cinderella`](https://sphire.mpg.de/wiki/doku.php?id=auto_2d_class_selection) - Cinderella: Deep learning based binary classification tool
* ### New clustering methods
1. [`γ-SUP`](http://sabid.stat.sinica.edu.tw/start) - [**γ-SUP: A clustering algorithm for cryo-electron microscopy images of asymmetric particles**](https://projecteuclid.org/euclid.aoas/1396966286)
2. [`Wasserstein-k-means`](https://github.com/4tywon/wasserstein-k-means) - [**Wasserstein K-Means for Clustering Tomographic Projections**](https://arxiv.org/abs/2010.09989)

## [3D Tomographic Reconstruction](https://www.sciencedirect.com/science/article/pii/S0076687910820014)
1. [**A Survey of the Use of Iterative Reconstruction Algorithms in Electron Microscopy**](https://www.ncbi.nlm.nih.gov/pmc/articles/PMC5623807/)
2. [`ASTRA`](https://github.com/astra-toolbox/astra-toolbox)

## Ab-initial model
* ### Class-averages based method
1. [`Simple`](https://simplecryoem.com/)-[**Single-particle cryo-EM-Improved Ab Initio 3D Reconstruction With SIMPLE/PRIME**](https://pubmed.ncbi.nlm.nih.gov/28795512/)
* ### Particles based method
1. `CryoSparc` - [**cryoSPARC: algorithms for rapid unsupervised cryo-EM structure determination**](https://www.nature.com/articles/nmeth.4169). [[Slides]](https://nramm.nysbc.org/wp-content/seminars/2017/slides/nysbc-nov2017-Brubaker.pdf)

## 3D refinement
* ### 3D Homogeneous Refinemnet
1. `Relion` - [**RELION: Implementation of a Bayesian approach to cryo-EM structure determination**](https://www.sciencedirect.com/science/article/pii/S1047847712002481). [[Video]](https://www.youtube.com/watch?time_continue=2&v=TfLFCeehfjM&feature=emb_title)
2. [`CryoSparc`](https://www.nature.com/articles/nmeth.4169) - Use Expectation-Maximization with branch and bound method for higher resolution
3. [`OPUS-SSRI`](https://github.com/alncat/cryoem) - [**Sparseness and Smoothness Regularized Imaging for improving the resolution of Cryo-EM single-particle reconstruction**](https://www.pnas.org/content/118/2/e2013756118)

* ### [New 3D reconstruction paradigm](https://www.sciencedirect.com/science/article/pii/S1047847722000909)
1. [`Pose Estimation with VAE-GAN`](https://github.com/compSPI/article-CVPR-workshop-2020)-[**Estimation of Orientation and Camera Parameters from Cryo-Electron Microscopy Images with Variational Autoencoders and Generative Adversarial**](https://openaccess.thecvf.com/content_CVPRW_2020/html/w57/Miolane_Estimation_of_Orientation_and_Camera_Parameters_From_Cryo-Electron_Microscopy_Images_CVPRW_2020_paper.html). [[Related work]](https://arxiv.org/abs/2107.02958)
2. [`CryoGAN`](https://github.com/harshit-gupta-cor/CryoGAN) - [**CryoGAN: A New Reconstruction Paradigm for Single-Particle Cryo-EM Via Deep Adversarial Learning**](https://ieeexplore.ieee.org/stamp/stamp.jsp?arnumber=9483649)
3. [`Orientation recovery with Siamese neural network`](https://github.com/JelenaBanjac/protein-reconstruction) - [**Learning to recover orientations from projections in single-particle cryo-EM**](https://arxiv.org/abs/2104.06237)
4. [`CryoAI`](https://github.com/compSPI/cryoAI) - [**Amortized Inference of Poses for Ab Initio Reconstruction of 3D Molecular Volumes from Real Cryo-EM Images**](https://arxiv.org/abs/2203.08138)

* ### 3D classification
1. [`Relion`](http://franklab.cpmc.columbia.edu/franklab/Learning_Materials/Meeting_7_Maximum_Likelihood/Papers/Scheres_2016_MIE_v579.pdf), `CryoSparc`  - Perturb the initial model and use projection matching with weighted assignment
2. [`LCTC`](https://github.com/ghl1995/LCTC) - [**An Efficient Method to Quantify Structural Distributions in Heterogeneous cryo-EM Datasets**](https://www.biorxiv.org/content/10.1101/2021.05.27.446075v1)

* ### 3D non-uniform Refinemnet
1. `CryoSparc` - [**Non-uniform refinement: Adaptive regularization improves single particle cryo-EM reconstruction**](https://www.nature.com/articles/s41592-020-00990-8)
2. [`SideSplitter`](https://github.com/StructuralBiology-ICLMedicine/SIDESPLITTER) - [**Mitigating Local Over-fitting During Single Particle Reconstruction with SIDESPLITTER**](https://www.sciencedirect.com/science/article/pii/S1047847720301180). [[Video]](https://www.youtube.com/watch?v=jTNH6Z0n254&list=PLFEB3YHuxu11Jp_pOCIEtXxSqozFHve0O&index=12)

## [3D variability analysis](http://scripts.iucr.org/cgi-bin/paper?S2053230X18015108)
1. [`sbackprop`](https://github.com/dkimanius/sbackprop) - [**Sparse Fourier Backpropagation in Cryo-EM Reconstruction**](https://openreview.net/pdf?id=51f5sPXJD_E)
1. [`Flexutils`](https://github.com/scipion-em/scipion-em-flexutils) - [**Estimating conformational landscapes from
Cryo-EM particles by 3D Zernike polynomials**](https://www.nature.com/articles/s41467-023-35791-y) ( Use in `Xmipp`)
1. `Diffusion Prior` - [**Latent Space Diffusion Models of Cryo-EM Structures**](https://arxiv.org/abs/2211.14169)
1. [`3DFlex`](https://cryosparc.com/3dflex) - [**3D Flexible Refinement: Structure and Motion of Flexible Proteins from Cryo-EM**](https://www.biorxiv.org/content/10.1101/2021.04.22.440893v1) ( Use in `CryoSparc`)
1. [`e2gmm`](https://blake.bcm.edu/emanwiki/EMAN2/e2gmm) - [**Deep learning-based mixed-dimensional Gaussian mixture model for characterizing variability in cryo-EM**](https://www.nature.com/articles/s41592-021-01220-5)  ( Use in `EMAN2`, more information see [here](https://arxiv.org/ftp/arxiv/papers/2211/2211.10518.pdf))
1. `3DVA` - [**3D Variability Analysis: Directly resolving continuous flexibility and discrete heterogeneity from single particle cryo-EM images**](https://www.sciencedirect.com/science/article/pii/S1047847721000071) ( Use in `CryoSparc`)
1. [`CryoDRGN2`](https://github.com/zhonge/cryodrgn/tree/vb/v2_quickstart) - [**CryoDRGN2: Ab initio neural reconstruction of 3D protein structures
from real cryo-EM images**](https://openaccess.thecvf.com/content/ICCV2021/papers/Zhong_CryoDRGN2_Ab_Initio_Neural_Reconstruction_of_3D_Protein_Structures_From_ICCV_2021_paper.pdf)
1. [`CryoDRGN`](https://github.com/zhonge/cryodrgn) - [**CryoDRGN: Reconstruction of heterogeneous structures from cryo-electron micrographs using neural networks**](https://www.nature.com/articles/s41592-020-01049-4?utm_source=other&utm_medium=other&utm_content=null&utm_campaign=JRCN_1_DD01_CN_NatureRJ_article_paid_XMOL) (For processing with large dataset see [here](https://zhonge.github.io/cryodrgn/index.html))
1. [`ManifoldEM`](https://github.com/evanseitz/ManifoldEM_Python) - [**Retrieving functional pathways of biomolecules from single-particle snapshots**](https://www.nature.com/articles/s41467-020-18403-x)
1. [`DMSA`](https://github.com/evanseitz/cryoEM_DMSA) - [**Recovery of conformational continuum from single-particle cryo-EM data: Optimization of ManifoldEM informed by ground-truth studies**](https://www.biorxiv.org/content/10.1101/2021.06.18.449029v1)
1. [`AlphaCryo4D`](https://github.com/AlphaCryo4D/AlphaCryo4D) - [**Visualizing Conformational Space of Functional Biomolecular Complexes by Deep Manifold Learning**](https://www.mdpi.com/1422-0067/23/16/8872)
1. [`BioEM`](https://github.com/bio-phys/BioEM) - [**A Bayesian approach to extracting free-energy profiles from cryo-electron microscopy experiments**](https://www.nature.com/articles/s41598-021-92621-1)
1. `VAE` - [**Inferring a Continuous Distribution of Atom Coordinates from Cryo-EM Images using VAEs**](https://arxiv.org/abs/2106.14108)
1. `cryoFIRE` - [**Amortized Inference for Heterogeneous Reconstruction in Cryo-EM**](https://arxiv.org/abs/2210.07387)
1. `Atomic VAE` - [**Heterogeneous reconstruction of deformable atomic models in Cryo-EM**](https://arxiv.org/abs/2209.15121)

* ### Latent space analysis
1. [`Polaris`](https://github.com/evanseitz/POLARIS) - [**POLARIS: Path of Least Action Analysis on Energy Landscapes**](https://pubs.acs.org/doi/full/10.1021/acs.jcim.9b01108)

* ### [Focus classifcation](https://www.frontiersin.org/articles/10.3389/fmolb.2019.00033/full)
1. [`localrec`](https://github.com/OPIC-Oxford/localrec) - [**Localized reconstruction of subunits from electron cryomicroscopy images of macromolecular complexes**](https://www.nature.com/articles/ncomms9843) (Use in `Scipion`)
2. `Multi-body refinement` - [**Characterisation of molecular motions in cryo-EM single-particle data by multi-body refinement in RELION**](https://elifesciences.org/articles/36861) (Use in `Relion`)

## [Postprocessing](https://www.sciencedirect.com/science/article/pii/S0079610720300699)
* ### Per-particle based motion and ctf refinement
1. `CTF refinement` - [`Relion3`](https://www.ncbi.nlm.nih.gov/pmc/articles/PMC6250425/)/[`CisTEM`](https://cistem.org/documentation#tab-1-12)/[`CryoSparc`](https://cryosparc.com/docs/tutorials/ctf-refinement), 3D Reference required
2. `Ewald sphere correction` - [**New tools for automated high-resolution cryo-EM structure determination in RELION-3**](https://elifesciences.org/articles/42166#bib45) (Use ` relion_reconstruct --reverse_curvature`)
3.  `High-order aberrations` - [**Estimation of high-order aberrations and anisotropic magnification from cryo-EM data sets in RELION-3.1**](https://journals.iucr.org/m/issues/2020/02/00/fq5009/index.html)
4. `Bayesian Polishing` - [**A Bayesian approach to beam-induced motion correction in cryo-EM single-particle analysis**](https://journals.iucr.org/m/issues/2019/01/00/fq5003/)  (Use in Relion, 3D Reference required)
5. [`M`](http://www.warpem.com/warp/?page_id=827) - [**Multi-particle cryo-EM refinement with M visualizes ribosome-antibiotic complex at 3.5 Å inside cells**](https://www.nature.com/articles/s41592-020-01054-7). [[Video]](https://www.youtube.com/watch?v=kiS-ELvQ1gc&list=PLFEB3YHuxu11Jp_pOCIEtXxSqozFHve0O&index=18)

* ### [Local resolution](https://www.sciencedirect.com/science/article/pii/S0959440X20301044) and directional resolution
1. [`Blocres`](https://lsbr.niams.nih.gov/bsoft/programs/blocres.html) - [**One number does not fit all: Mapping local variations in resolution in cryo-EM reconstructions**](https://www.ncbi.nlm.nih.gov/pmc/articles/PMC3837392/)
1. [`ResMap`](https://sourceforge.net/projects/resmap-latest/) - [**Quantifying the Local Resolution of Cryo-EM Density Maps**](https://www.nature.com/articles/nmeth.2727)
2. `MonoRes` - [**MonoRes: Automatic and Accurate Estimation of Local Resolution for Electron Microscopy Maps**](https://pubmed.ncbi.nlm.nih.gov/29395788/) (Use in `Scipion`)
3. [`3DFSC`](https://github.com/nysbc/Anisotropy) - [**Addressing preferred specimen orientation in single-particle cryo-EM through tilting**](https://www.nature.com/articles/nmeth.4347)
4. `MonoDir` - [**Measuring local-directional resolution and local anisotropy in cryo-EM maps**](https://www.nature.com/articles/s41467-019-13742-w) (Use in `Scipion`)

* ### [Sharpening](https://www.sciencedirect.com/science/article/pii/S0022283603010222) and local filtering to enhance [interpretability](https://www.sciencedirect.com/science/article/pii/S0079610720300699)
1. [`EMReady`](http://huanglab.phys.hust.edu.cn/EMReady/)
2. [`Confidence Maps`](https://git.embl.de/mbeckers/FDRthresholding/-/tree/master) - [**Thresholding of cryo-EM density maps by false discovery rate control**](https://journals.iucr.org/m/issues/2019/01/00/pw5002/index.html) (Use in [`ccpem`](https://www.ccpem.ac.uk/download.php))
3. [`DeepEnhancer`](https://github.com/rsanchezgarc/deepEMhancer)-[**DeepEMhancer: a deep learning solution for cryo-EM volume post-processing**](https://www.biorxiv.org/content/10.1101/2020.06.12.148296v3)
4. [`LocScale`](https://git.embl.de/jakobi/LocScale)[**Model-based local density sharpening of cryo-EM maps**](https://elifesciences.org/articles/27131) (Use in [`ccpem`](https://www.ccpem.ac.uk/download.php))
5. [`LocalDeblur`](https://github.com/I2PC/scipion-em-xmipp/wiki/XmippProtLocSharp) - [**Automatic local resolution-based sharpening of cryo-EM maps**](https://academic.oup.com/bioinformatics/article/36/3/765/5554698)
6. [`EM-GAN`](https://kiharalab.org/emsuites/emgan.php) - **Improved Protein Structure Modeling Using Enhanced Cryo-EM Maps With 3D Deep Generative Networks**

* ### Denoising 3D volume
1. [`Warp`](http://www.warpem.com/warp/?page_id=389) - Based on Noise2Noise
2. `Topaz` - Contains 3D denoise functionality
3. [`Relion`](https://github.com/3dem/externprior) - [**Exploiting prior knowledge about biological macromolecules in cryo-EM structure determination**](https://journals.iucr.org/m/issues/2021/01/00/fq5015/)

* ### Map to model
1. [`ModelAngelo`](https://github.com/3dem/model-angelo) - [**ModelAngelo: Automated Model Building in Cryo-EM Maps**](https://arxiv.org/abs/2210.00006)

## Visualization
1. [`Chimera`](https://www.cgl.ucsf.edu/chimera/) - [**UCSF Chimera--a Visualization System for Exploratory Research and Analysis**](https://pubmed.ncbi.nlm.nih.gov/15264254/)
2. [`ChimeraX`](http://www.rbvi.ucsf.edu/chimerax/) - [**UCSF ChimeraX: Meeting modern challenges in visualization and analysis**](https://onlinelibrary.wiley.com/doi/full/10.1002/pro.3235)

# Conventions
## 3DEM Convention
* [3DEM Conventions](https://github.com/azazellochg/3DEM-conventions)
* [Resolution Measure](https://www.sciencedirect.com/science/article/pii/S0079610716300037)
* [Symmetry](https://scipion-em.github.io/docs/docs/developer/symmetries)

## Image contrast
`White on Black` - Relion, Xmipp, EMAN2, CryoSparc (Cryo-EM data is typically recorded as Black on White but will invert during processing)

`Black on White` - Frealign (CisTEM)

## Pixel size
For the format defined in RELION, the actual pixel size is calculated as `rlnDetectorPixelSize * 10000 / rlnMagnification`.
In Relion 3.1 it has been replaced with `rlnImagePixelSize`.

## FSC calculation
Using [e2proc3d](http://sphire.mpg.de/wiki/doku.php?id=pipeline:utilities:e2proc3d)

## Mask generation in [CryoSparc](https://discuss.cryosparc.com/t/tight-corrected-and-loose-gsfsc-curves/201/5)

### No Mask:
This is the raw FSC calculated between two independent half-maps reconstructed from the data. There is no masking applied, so both the structure and solvent are included in this FSC.

### Spherical:
This is the FSC calculated after applying a soft spherical mask to both half maps. The outer radius of the soft sphere is equal to half the volume box-size (i.e. the sphere extends to the faces of the box in all directions). The inner radius is 85 percent of the outer radius. Between inner and outer radii, a soft cosine edge transitions from a mask value of one to a value of zero.

### Loose:
This is the FSC calculated after applying a soft solvent mask to both half maps. The loose mask is calculated as follows. First, the density map is thresholded at 50% of the maximum density value. The resulting volume is dilated to create a soft mask. Voxels in the mask that are within 25 angstroms of the thresholded region receive a mask value of 1.0. Voxels between 25 and 40 angstroms fall off with a soft cosine edge, and voxels outside 40 angstroms receive a value of 0.0.

### Tight:
This is the same as the loose mask, except the dilation distances are 6 angstroms for the value 1.0 distance and 12 angstroms for the value 0.0 distance.

### Corrected:
This is the FSC curve calculated using the tight mask with correction by noise substitution [1]. The two half maps have their phases randomized beyond a certain resolution, then the tight mask is applied to both, and an FSC is calculated. This FSC is used along with the original FSC before phase randomization to compute the corrected FSC as in [1]. This accounts for correlation effects induced by masking. The resolution at which phase randomization begins is the resolution at which the no-mask FSC drops below the FSC = 0.143 criteria.

> Chen, S. et al. High-resolution noise substitution to measure overfitting and validate resolution in 3D structure determination by single particle electron cryomicroscopy. Ultramicroscopy 135, 24–35 (2013).

## Sharpen and filtering in [CryoSparc](https://discuss.cryosparc.com/t/how-is-the-non-uniform-refinement-map-filtered-map-created/4388)
The `map_filtered` output in non-uniform refinement is generated as follows, after refinement has converged:

1. both raw, unfiltered halfmaps are averaged together
2. the raw map is filtered using the Gold-Standard FSC curve
3. the filtered map is sharpened using the Guinier-plot estimated b-factor
4. the sharpened map is outputted as `_map_sharp_local.mrc` (confusing filename… sorry)
5. the sharpened map is locally filtered using a local resolution estimate computed from the half-maps and the locally filtered map is outputted as `_map_filtered.mrc`

# Tips
## Format Conversion
1. Using [PyRelion](https://github.com/OPIC-Oxford/localrec/wiki/Tips-and-tricks)
Operate Star files

1. Using [PyEM](https://github.com/asarnow/pyem)
Convert metadata

1. Using [EMAN2](https://github.com/cryoem/eman2/blob/master/examples/e2convertrelion.py)
Convert from relion to eman2

1. Using [EMAN2](https://blake.bcm.edu/emanwiki/EMAN2)
Convert binary data

## Parse `Star file`
1. Using [`starfile`](https://github.com/teamtomo/starfile) or [`starparser`](https://github.com/sami-chaaban/starparser)
2. Using Custom function like [here](Parse_Star.ipynb)

## Parse `CryoSparc file`
1. Using [`PyEM`](https://github.com/asarnow/pyem)
2. Using Custom function like [here](Parse_Cryosparc.ipynb). See full tutorial [here](https://guide.cryosparc.com/processing-data/manipulating-.cs-files-created-by-cryosparc)
3. [`cryosparc-tools`](https://github.com/cryoem-uoft/cryosparc-tools)

## Performing Focus Classification
### With Chimera
### [CisTEM](https://cistem.org/frequently-asked-questions#tab-1-4)
A focus mask is defined as a sphere specified by radius and x,y,z coordinates of the sphere center, all given in Å. In order to find x,y,z, in CisTEM
1. Open the 3D map in Chimera. The map should be assigned model #0.
1. Display the command line by selecting Tools > General Controls > Command Line.
1. Execute vop threshold #0 maximum -1000 setMaximum 1.0 (this assumes that the minimum voxel value in the 3D map is larger than 1000). This will create a new model #1.
1. Execute shape sphere radius 30 color red mesh true on the Command Line. This will generate a sphere with a 30 Å radius, numbered model #2 (change the radius as needed).
1. Open Tools > General Controls > Model Panel and deactivate models #0 and #1. Then move the sphere (model #2) into model #1 until it is fully contained in that volume.
1. Execute mask #1 #2 on the Command Line. This will generate a new model #3 of a solid sphere with a 30 Å radius.
1. Close models #1 and #2 in Model Panel and make sure model #0 is visible but deactive.
1. Select model #3 (the solid sphere) in Volume Viewer and change its color to red (or something easily distiguishable from the original 3D map) and make it transparent.
1. Move model #3 to the desired location of the focus mask in the 3D map.
1. Execute vop resample #3 onGrid #0 on the Command Line. This will generate a copy of the solid sphere as a new model #1, now resampled in the same coordinate system as the original 3D map.
1. Execute measure center #1 on the Command Line. This will display the x,y,z coordinates of the mask in pixel coordinates below the command line. These coordinates have to be converted to Å by multiplying them with the pixel size of the 3D map before they can be used in cisTEM's Manual Refine panel.

### [Local classification in CryoSparc](https://cryosparc.com/blog/local-refinement-snRNP-case-study/)

## Display images
### Micrographs
In Matlab use `grayImage = uint8(255 * mat2gray(originalImage)); imshow(grayImage);`

## Using Chimera
### Create synthetic map from PDB
1. Open chimera and import that protein: File -> Fetch By ID… -> Select PDB and type the PDB ID of the protein -> Fetch. Now you will see the 3D structure of the protein.
2. Create synthetic map from PDB model: Modify the density of the protein: Tools -> General Controls -> Command Line -> Type: molmap #0 5 (molmap = is command that generates a density map from the specified atom; # is atom specification, i.e. number assigned to the model by default; 5 = resolution).
3. Store reference map as MRC: Save to a file: Tools -> Volume Data -> Volume Viewer -> File -> Save map as… -> Give it the protein_PDB_ID.mrc.

### [Dispaly local resolution](https://discuss.cryosparc.com/t/view-coloured-3d-map-in-chimera-using-output-from-local-resolution/1003)
### [Flip handedness of model and map](https://github.com/asarnow/pyem/wiki/Z-Flip-for-Models-and-Maps)
### Slice view
`volume #0 planes z,220 step 1 level -1 style surface`
### [Visualizing dynamic movies](https://guide.cryosparc.com/processing-data/tutorials-and-case-studies/tutorial-3d-variability-analysis-part-one#visualization)
#### [Using Script](https://github.com/zhonge/cryodrgn/discussions/79)

# Contributing
Your contributions are highly appreciated! Please take a look at the [contribution guidelines](https://github.com/phonchi/Computational-CryoEM/blob/master/CONTRIBUTING.md) first.


text:
include distribute_setup.py
include REQUIREMENTS
include REQUIREMENTS.*
include README.md
include UNLICENSE

repo_name: colinhowe/pyFaceGraph

repo_language: Python


text:
#!python
"""Bootstrap distribute installation

If you want to use setuptools in your package's setup.py, just include this
file in the same directory with it, and add this to the top of your setup.py::

    from distribute_setup import use_setuptools
    use_setuptools()

If you want to require a specific version of setuptools, set a download
mirror, or use an alternate download directory, you can do so by supplying
the appropriate options to ``use_setuptools()``.

This file can also be run as a script to install or upgrade setuptools.
"""
import os
import sys
import time
import fnmatch
import tempfile
import tarfile
from distutils import log

try:
    from site import USER_SITE
except ImportError:
    USER_SITE = None

try:
    import subprocess

    def _python_cmd(*args):
        args = (sys.executable,) + args
        return subprocess.call(args) == 0

except ImportError:
    # will be used for python 2.3
    def _python_cmd(*args):
        args = (sys.executable,) + args
        # quoting arguments if windows
        if sys.platform == 'win32':
            def quote(arg):
                if ' ' in arg:
                    return '"%s"' % arg
                return arg
            args = [quote(arg) for arg in args]
        return os.spawnl(os.P_WAIT, sys.executable, *args) == 0

DEFAULT_VERSION = "0.6.10"
DEFAULT_URL = "http://pypi.python.org/packages/source/d/distribute/"
SETUPTOOLS_FAKED_VERSION = "0.6c11"

SETUPTOOLS_PKG_INFO = """\
Metadata-Version: 1.0
Name: setuptools
Version: %s
Summary: xxxx
Home-page: xxx
Author: xxx
Author-email: xxx
License: xxx
Description: xxx
""" % SETUPTOOLS_FAKED_VERSION

def _install(tarball):
    # extracting the tarball
    tmpdir = tempfile.mkdtemp()
    log.warn('Extracting in %s', tmpdir)
    old_wd = os.getcwd()
    try:
        os.chdir(tmpdir)
        tar = tarfile.open(tarball)
        _extractall(tar)
        tar.close()

        # going in the directory
        subdir = os.path.join(tmpdir, os.listdir(tmpdir)[0])
        os.chdir(subdir)
        log.warn('Now working in %s', subdir)

        # installing
        log.warn('Installing Distribute')
        if not _python_cmd('setup.py', 'install'):
            log.warn('Something went wrong during the installation.')
            log.warn('See the error message above.')
    finally:
        os.chdir(old_wd)

def _build_egg(egg, tarball, to_dir):
    # extracting the tarball
    tmpdir = tempfile.mkdtemp()
    log.warn('Extracting in %s', tmpdir)
    old_wd = os.getcwd()
    try:
        os.chdir(tmpdir)
        tar = tarfile.open(tarball)
        _extractall(tar)
        tar.close()

        # going in the directory
        subdir = os.path.join(tmpdir, os.listdir(tmpdir)[0])
        os.chdir(subdir)
        log.warn('Now working in %s', subdir)

        # building an egg
        log.warn('Building a Distribute egg in %s', to_dir)
        _python_cmd('setup.py', '-q', 'bdist_egg', '--dist-dir', to_dir)

    finally:
        os.chdir(old_wd)
    # returning the result
    log.warn(egg)
    if not os.path.exists(egg):
        raise IOError('Could not build the egg.')

def _do_download(version, download_base, to_dir, download_delay):
    egg = os.path.join(to_dir, 'distribute-%s-py%d.%d.egg'
                       % (version, sys.version_info[0], sys.version_info[1]))
    if not os.path.exists(egg):
        tarball = download_setuptools(version, download_base,
                                      to_dir, download_delay)
        _build_egg(egg, tarball, to_dir)
    sys.path.insert(0, egg)
    import setuptools
    setuptools.bootstrap_install_from = egg

def use_setuptools(version=DEFAULT_VERSION, download_base=DEFAULT_URL,
                   to_dir=os.curdir, download_delay=15, no_fake=True):
    # making sure we use the absolute path
    to_dir = os.path.abspath(to_dir)
    was_imported = 'pkg_resources' in sys.modules or \
        'setuptools' in sys.modules
    try:
        try:
            import pkg_resources
            if not hasattr(pkg_resources, '_distribute'):
                if not no_fake:
                    _fake_setuptools()
                raise ImportError
        except ImportError:
            return _do_download(version, download_base, to_dir, download_delay)
        try:
            pkg_resources.require("distribute>="+version)
            return
        except pkg_resources.VersionConflict:
            e = sys.exc_info()[1]
            if was_imported:
                sys.stderr.write(
                "The required version of distribute (>=%s) is not available,\n"
                "and can't be installed while this script is running. Please\n"
                "install a more recent version first, using\n"
                "'easy_install -U distribute'."
                "\n\n(Currently using %r)\n" % (version, e.args[0]))
                sys.exit(2)
            else:
                del pkg_resources, sys.modules['pkg_resources']    # reload ok
                return _do_download(version, download_base, to_dir,
                                    download_delay)
        except pkg_resources.DistributionNotFound:
            return _do_download(version, download_base, to_dir,
                                download_delay)
    finally:
        if not no_fake:
            _create_fake_setuptools_pkg_info(to_dir)

def download_setuptools(version=DEFAULT_VERSION, download_base=DEFAULT_URL,
                        to_dir=os.curdir, delay=15):
    """Download distribute from a specified location and return its filename

    `version` should be a valid distribute version number that is available
    as an egg for download under the `download_base` URL (which should end
    with a '/'). `to_dir` is the directory where the egg will be downloaded.
    `delay` is the number of seconds to pause before an actual download
    attempt.
    """
    # making sure we use the absolute path
    to_dir = os.path.abspath(to_dir)
    try:
        from urllib.request import urlopen
    except ImportError:
        from urllib2 import urlopen
    tgz_name = "distribute-%s.tar.gz" % version
    url = download_base + tgz_name
    saveto = os.path.join(to_dir, tgz_name)
    src = dst = None
    if not os.path.exists(saveto):  # Avoid repeated downloads
        try:
            log.warn("Downloading %s", url)
            src = urlopen(url)
            # Read/write all in one block, so we don't create a corrupt file
            # if the download is interrupted.
            data = src.read()
            dst = open(saveto, "wb")
            dst.write(data)
        finally:
            if src:
                src.close()
            if dst:
                dst.close()
    return os.path.realpath(saveto)

def _patch_file(path, content):
    """Will backup the file then patch it"""
    existing_content = open(path).read()
    if existing_content == content:
        # already patched
        log.warn('Already patched.')
        return False
    log.warn('Patching...')
    _rename_path(path)
    f = open(path, 'w')
    try:
        f.write(content)
    finally:
        f.close()
    return True

def _same_content(path, content):
    return open(path).read() == content

def _no_sandbox(function):
    def __no_sandbox(*args, **kw):
        try:
            from setuptools.sandbox import DirectorySandbox
            def violation(*args):
                pass
            DirectorySandbox._old = DirectorySandbox._violation
            DirectorySandbox._violation = violation
            patched = True
        except ImportError:
            patched = False

        try:
            return function(*args, **kw)
        finally:
            if patched:
                DirectorySandbox._violation = DirectorySandbox._old
                del DirectorySandbox._old

    return __no_sandbox

@_no_sandbox
def _rename_path(path):
    new_name = path + '.OLD.%s' % time.time()
    log.warn('Renaming %s into %s', path, new_name)
    os.rename(path, new_name)
    return new_name

def _remove_flat_installation(placeholder):
    if not os.path.isdir(placeholder):
        log.warn('Unkown installation at %s', placeholder)
        return False
    found = False
    for file in os.listdir(placeholder):
        if fnmatch.fnmatch(file, 'setuptools*.egg-info'):
            found = True
            break
    if not found:
        log.warn('Could not locate setuptools*.egg-info')
        return

    log.warn('Removing elements out of the way...')
    pkg_info = os.path.join(placeholder, file)
    if os.path.isdir(pkg_info):
        patched = _patch_egg_dir(pkg_info)
    else:
        patched = _patch_file(pkg_info, SETUPTOOLS_PKG_INFO)

    if not patched:
        log.warn('%s already patched.', pkg_info)
        return False
    # now let's move the files out of the way
    for element in ('setuptools', 'pkg_resources.py', 'site.py'):
        element = os.path.join(placeholder, element)
        if os.path.exists(element):
            _rename_path(element)
        else:
            log.warn('Could not find the %s element of the '
                     'Setuptools distribution', element)
    return True

def _after_install(dist):
    log.warn('After install bootstrap.')
    placeholder = dist.get_command_obj('install').install_purelib
    _create_fake_setuptools_pkg_info(placeholder)

@_no_sandbox
def _create_fake_setuptools_pkg_info(placeholder):
    if not placeholder or not os.path.exists(placeholder):
        log.warn('Could not find the install location')
        return
    pyver = '%s.%s' % (sys.version_info[0], sys.version_info[1])
    setuptools_file = 'setuptools-%s-py%s.egg-info' % \
            (SETUPTOOLS_FAKED_VERSION, pyver)
    pkg_info = os.path.join(placeholder, setuptools_file)
    if os.path.exists(pkg_info):
        log.warn('%s already exists', pkg_info)
        return

    log.warn('Creating %s', pkg_info)
    f = open(pkg_info, 'w')
    try:
        f.write(SETUPTOOLS_PKG_INFO)
    finally:
        f.close()

    pth_file = os.path.join(placeholder, 'setuptools.pth')
    log.warn('Creating %s', pth_file)
    f = open(pth_file, 'w')
    try:
        f.write(os.path.join(os.curdir, setuptools_file))
    finally:
        f.close()

def _patch_egg_dir(path):
    # let's check if it's already patched
    pkg_info = os.path.join(path, 'EGG-INFO', 'PKG-INFO')
    if os.path.exists(pkg_info):
        if _same_content(pkg_info, SETUPTOOLS_PKG_INFO):
            log.warn('%s already patched.', pkg_info)
            return False
    _rename_path(path)
    os.mkdir(path)
    os.mkdir(os.path.join(path, 'EGG-INFO'))
    pkg_info = os.path.join(path, 'EGG-INFO', 'PKG-INFO')
    f = open(pkg_info, 'w')
    try:
        f.write(SETUPTOOLS_PKG_INFO)
    finally:
        f.close()
    return True

def _before_install():
    log.warn('Before install bootstrap.')
    _fake_setuptools()

def _under_prefix(location):
    if 'install' not in sys.argv:
        return True
    args = sys.argv[sys.argv.index('install')+1:]
    for index, arg in enumerate(args):
        for option in ('--root', '--prefix'):
            if arg.startswith('%s=' % option):
                top_dir = arg.split('root=')[-1]
                return location.startswith(top_dir)
            elif arg == option:
                if len(args) > index:
                    top_dir = args[index+1]
                    return location.startswith(top_dir)
            elif option == '--user' and USER_SITE is not None:
                return location.startswith(USER_SITE)
    return True

def _fake_setuptools():
    log.warn('Scanning installed packages')
    try:
        import pkg_resources
    except ImportError:
        # we're cool
        log.warn('Setuptools or Distribute does not seem to be installed.')
        return
    ws = pkg_resources.working_set
    try:
        setuptools_dist = ws.find(pkg_resources.Requirement.parse('setuptools',
                                  replacement=False))
    except TypeError:
        # old distribute API
        setuptools_dist = ws.find(pkg_resources.Requirement.parse('setuptools'))

    if setuptools_dist is None:
        log.warn('No setuptools distribution found')
        return
    # detecting if it was already faked
    setuptools_location = setuptools_dist.location
    log.warn('Setuptools installation detected at %s', setuptools_location)

    # if --root or --preix was provided, and if
    # setuptools is not located in them, we don't patch it
    if not _under_prefix(setuptools_location):
        log.warn('Not patching, --root or --prefix is installing Distribute'
                 ' in another location')
        return

    # let's see if its an egg
    if not setuptools_location.endswith('.egg'):
        log.warn('Non-egg installation')
        res = _remove_flat_installation(setuptools_location)
        if not res:
            return
    else:
        log.warn('Egg installation')
        pkg_info = os.path.join(setuptools_location, 'EGG-INFO', 'PKG-INFO')
        if (os.path.exists(pkg_info) and
            _same_content(pkg_info, SETUPTOOLS_PKG_INFO)):
            log.warn('Already patched.')
            return
        log.warn('Patching...')
        # let's create a fake egg replacing setuptools one
        res = _patch_egg_dir(setuptools_location)
        if not res:
            return
    log.warn('Patched done.')
    _relaunch()

def _relaunch():
    log.warn('Relaunching...')
    # we have to relaunch the process
    args = [sys.executable] + sys.argv
    sys.exit(subprocess.call(args))

def _extractall(self, path=".", members=None):
    """Extract all members from the archive to the current working
       directory and set owner, modification time and permissions on
       directories afterwards. `path' specifies a different directory
       to extract to. `members' is optional and must be a subset of the
       list returned by getmembers().
    """
    import copy
    import operator
    from tarfile import ExtractError
    directories = []

    if members is None:
        members = self

    for tarinfo in members:
        if tarinfo.isdir():
            # Extract directories with a safe mode.
            directories.append(tarinfo)
            tarinfo = copy.copy(tarinfo)
            tarinfo.mode = 448 # decimal for oct 0700
        self.extract(tarinfo, path)

    # Reverse sort directories.
    if sys.version_info < (2, 4):
        def sorter(dir1, dir2):
            return cmp(dir1.name, dir2.name)
        directories.sort(sorter)
        directories.reverse()
    else:
        directories.sort(key=operator.attrgetter('name'), reverse=True)

    # Set correct owner, mtime and filemode on directories.
    for tarinfo in directories:
        dirpath = os.path.join(path, tarinfo.name)
        try:
            self.chown(tarinfo, dirpath)
            self.utime(tarinfo, dirpath)
            self.chmod(tarinfo, dirpath)
        except ExtractError:
            e = sys.exc_info()[1]
            if self.errorlevel > 1:
                raise
            else:
                self._dbg(1, "tarfile: %s" % e)

def main(argv, version=DEFAULT_VERSION):
    """Install or upgrade setuptools and EasyInstall"""
    tarball = download_setuptools()
    _install(tarball)

if __name__ == '__main__':
    main(sys.argv[1:])


text:
bunch>=1.0.1
ujson>=1.19


text:
#!/usr/bin/env python
# -*- coding: utf-8 -*-

import glob
import os
import re

from distribute_setup import use_setuptools; use_setuptools()
from setuptools import setup, find_packages

rel_file = lambda *args: os.path.join(os.path.dirname(os.path.abspath(__file__)), *args)
cleanup = lambda lines: filter(None, map(lambda s: s.strip(), lines))

def read_from(filename):
    fp = open(filename)
    try:
        return fp.read()
    finally:
        fp.close()

def get_version():
    data = read_from(rel_file('src', 'facegraph', '__init__.py'))
    return re.search(r"__version__ = '([^']+)'", data).group(1)

def get_requirements():
    return cleanup(read_from(rel_file('REQUIREMENTS')).splitlines())

def get_extra_requirements():
    extras_require = {}
    for req_filename in glob.glob(rel_file('REQUIREMENTS.*')):
        group = os.path.basename(req_filename).split('.', 1)[1]
        extras_require[group] = cleanup(read_from(req_filename).splitlines())
    return extras_require

setup(
    name             = 'pyfacegraph',
    version          = get_version(),
    author           = "iPlatform Ltd",
    author_email     = "opensource@theiplatform.com",
    url              = 'http://github.com/iplatform/pyFaceGraph',
    description      = "A client library for the Facebook Graph API.",
    packages         = find_packages(where='src'),
    package_dir      = {'': 'src'},
    install_requires = get_requirements(),
    extras_require   = get_extra_requirements(),
)


text:
# pyFaceGraph

pyFaceGraph is a Python client library for the [Facebook Graph API][api]. It is
being developed and maintained by [iPlatform][].

  [api]: http://developers.facebook.com/docs/api
  [iplatform]: http://theiplatform.com/

## Installation

Via `pip` or `easy_install`:

    pip install pyfacegraph
    easy_install pyfacegraph

You can install an 'edge' version via git:

    git clone 'git://github.com/iplatform/pyFaceGraph.git'
    cd pyFaceGraph
    python setup.py install

## Graph API Tutorial

To begin using the API, create a new `Graph` with an access token:

    >>> from facegraph import Graph
    >>> g = Graph(ACCESS_TOKEN) #  Access token is optional.
    >>> g
    <Graph('https://graph.facebook.com/') at 0x...>

### Addressing Nodes

Each `Graph` contains an access token and a URL. The graph you just created
will have a URL of 'https://graph.facebook.com/' by default (this is defined as
the class attribute `Graph.API_ROOT`).

    >>> print g.url
    https://graph.facebook.com/
    >>> unicode(g.url)
    u'https://graph.facebook.com/'

To address child nodes, `Graph` supports dynamic attribute and item lookups:

    >>> g.me
    <Graph('https://graph.facebook.com/me') at 0x...>
    >>> g.me.home
    <Graph('https://graph.facebook.com/me/home') at 0x...>
    >>> g['me']['home']
    <Graph('https://graph.facebook.com/me/home') at 0x...>
    >>> g[123456789]
    <Graph('https://graph.facebook.com/123456789') at 0x...>

Note that a `Graph` instance is rarely modified; these methods all return copies
of the original object. In addition, the API is lazy: HTTP requests will
never be made unless you explicitly make them.

### Retrieving Nodes

You can fetch data by calling a `Graph` instance:

    >>> about_me = g.me()
    >>> about_me
    Node({'about': '...', 'id': '1503223370'})

This returns a `Node` object, which contains the retrieved data. `Node` is
a subclass of `bunch.Bunch` [[pypi](http://pypi.python.org/pypi/bunch)], so you
can access keys using either attribute or item syntax:

    >>> about_me.id
    '1503223370'
    >>> about_me.first_name
    'Zachary'
    >>> about_me.hometown.name
    'London, United Kingdom'
    >>> about_me['hometown']['name']
    'London, United Kingdom'

Accessing non-existent attributes or items will return a `Graph` instance
corresponding to a child node. This `Graph` can then be called normally, to
retrieve the child node it represents:

    >>> 'home' in about_me  # Not present in the data itself
    False
    >>> about_me.home
    <Graph('https://graph.facebook.com/me/home') at 0x...>
    >>> about_me.home()
    Node({'data': [...]})

### Creating, Updating and Deleting Nodes

With the Graph API, node manipulation is done via HTTP POST requests.
`Graph.post()` will POST to the current URL, with varying semantics for each
endpoint:

    >>> post = g.me.feed.post(message="Test.")  # Status update
    >>> post
    Node({'id': '...'})
    >>> g[post.id].comments.post(message="A comment.") # Comment on the post
    Node({'id': '...'})
    >>> g[post.id].likes.post()  # Like the post
    True

    >>> event = g[121481007877204]()
    >>> event.name
    'Facebook Developer Garage London May 2010'
    >>> event.rsvp_status is None
    True
    >>> event.attending.post()  # Attend the given event
    True

Any keyword arguments passed to `post()` will be added as form data. Consult the
[Facebook API docs][api-ref] for a complete reference on URLs and options.

  [api-ref]: http://developers.facebook.com/docs/reference/api/

Nodes can be deleted by adding `?method=delete` to the URL; the `delete()`
method is a helpful shortcut:

    >>> g[post.id].delete()
    True

## Canvas Applications

Facebook is rolling out support for OAuth 2.0 in canvas applications with the
new `signed_request` parameter; you can read about it [here][signed_request].

  [signed_request]: http://developers.facebook.com/docs/authentication/canvas

pyFaceGraph comes with a simple function for verifying and decoding this
parameter:

    >>> from facegraph import decode_signed_request, InvalidSignature

    >>> decode_signed_request(APP_SECRET, GET['signed_request'])
    {'0': 'payload', 'algorithm': 'HMAC-SHA256'}

    >>> decode_signed_request('wrong-secret', GET['signed_request'])
    Traceback (most recent call last):
        ...
    InvalidSignature

### OAuth

pyFaceGraph defines two abstract class-based views (using [django-clsview][])
and a mixin; these help your app obtain OAuth 2.0 access tokens to access the
Graph API on behalf of a Facebook user.

  [django-clsview]: http://github.com/zacharyvoase/django-clsview

*   `FacebookOAuthView`: Defines some common methods useful to all
    Facebook-related CBVs. Should be used as an inherited ‘mixin’, i.e.:

        class Callback(CallbackView, FacebookOAuthView):
            pass

    Examples of methods implemented here are `redirect_uri()` (which must be the
    same for both the authorize and callback views), `client_id()` and
    `client_secret()` and `fetch_url()`

*   `AuthorizeView`: Acts solely to redirect users to Facebook for
    authorization. You can override `authorize_url()` to change the URL which
    the user is redirected to. `scope()` and `display()` are two shortcuts for
    the Facebook-specific parameters, which you might need to alter for your
    requirements.

*   `CallbackView`: The view the user is redirected back to from Facebook, upon
    successful authorization. It’s up to you to write the body of this (by
    defining `__call__()`), but you can call `get_access_token()` to fetch the
    access token. Suggestions include saving the access token to the database,
    or storing it in `request.session`.

See the source of [facegraph.django.views][fdv] for the ultimate reference guide
to these views.

  [fdv]: http://github.com/iplatform/pyFaceGraph/blob/master/src/facegraph/django/views.py

#### Example

Take a look at the example project in `test/graphdevtools`, which shows how to
construct an app that interfaces with Facebook (and provides a dashboard for
generating access tokens with varying permissions). You’ll need to create
`local_settings.py` in that directory and define `SECRET_KEY`,
`FACEBOOK_CLIENT_ID` and `FACEBOOK_CLIENT_SECRET` (for the latter two you need
to register a Facebook application).

### Middleware

The provided `FacebookGraphMiddleware` will attach a `facegraph.Graph` instance
to each request, accessible as `request.graph`. You will need to subclass this
middleware to define your own method of fetching the access token. For example:

    ## myapp/middleware.py:

    from facegraph.django.middleware import FacebookGraphMiddleware

    class GraphMiddleware(FacebookGraphMiddleware):
        def access_token(self, request):
            return request.session.get('access_token')

    ## settings.py:

    MIDDLEWARE_CLASSES = (
      # ...
      'myapp.middleware.GraphMiddleware',
      # ...
    )

Note that this will still attach a `Graph` even if the access token is `None`.
To check for authentication, just use `if request.graph.access_token:` in your
view code.

There is also a piece of middleware for canvas applications:
`FacebookCanvasMiddleware` will verify and decode Facebook’s new
`signed_request` parameter for each request, and attach the value as
`request.fbrequest`. By default, `settings.FACEBOOK_APP_SECRET` will be used as
the application secret for signature verification, but you can subclass the
middleware and override the `app_secret()` method to modify this behavior:

    ## myapp/middleware.py:

    from django.conf import settings
    from facegraph.django.middleware import FacebookCanvasMiddleware

    class CanvasMiddleware(FacebookCanvasMiddleware):
        def app_secret(self, request):
            return settings.MY_FACEBOOK_APP_SECRET

    ## settings.py:

    MIDDLEWARE_CLASSES = (
      # ...
      'myapp.middleware.CanvasMiddleware',
      # ...
    )

Finally, `FacebookCanvasGraphMiddleware` is a subclass of
`FacebookGraphMiddleware` which will use the signed request as the source of the
access token; no subclassing is necessary, just place it *after*
`CanvasMiddleware` in `MIDDLEWARE_CLASSES`:

    MIDDLEWARE_CLASSES = (
      # ...
      'myapp.middleware.CanvasMiddleware',
      'facegraph.django.middleware.FacebookCanvasGraphMiddleware',
      # ...
    )

## (Un)license

This is free and unencumbered software released into the public domain.

Anyone is free to copy, modify, publish, use, compile, sell, or distribute this
software, either in source code form or as a compiled binary, for any purpose,
commercial or non-commercial, and by any means.

In jurisdictions that recognize copyright laws, the author or authors of this
software dedicate any and all copyright interest in the software to the public
domain. We make this dedication for the benefit of the public at large and to
the detriment of our heirs and successors. We intend this dedication to be an
overt act of relinquishment in perpetuity of all present and future rights to
this software under copyright law.

THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS
FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS BE
LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF
CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE
SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.

For more information, please refer to <http://unlicense.org/>


text:
import urllib
import urlparse

def get_path(url):
    scheme, host, path, query, fragment = urlparse.urlsplit(url)
    return path

def get_host(url):
    scheme, host, path, query, fragment = urlparse.urlsplit(url)
    return host

def add_path(url, new_path):
    """Given a url and path, return a new url that combines
    the two.
    """
    scheme, host, path, query, fragment = urlparse.urlsplit(url)
    new_path = new_path.lstrip('/')
    if path.endswith('/'):
        path += new_path
    else:
        path += '/' + new_path
    return urlparse.urlunsplit([scheme, host, path, query, fragment])

def _query_param(key, value):
    """ensure that a query parameter's value is a string
    of bytes in UTF-8 encoding.
    """
    if isinstance(value, unicode):
        pass
    elif isinstance(value, str):
        value = value.decode('utf-8')
    else:
        value = unicode(value)
    return key, value.encode('utf-8')

def _make_query_tuples(params):
    if hasattr(params, 'items'):
        return [_query_param(*param) for param in params.items()]
    else:
        return [_query_param(*params)]

def add_query_params(url, params):
    """use the _update_query_params function to set a new query
    string for the url based on params.
    """
    return update_query_params(url, params, update=False)

def update_query_params(url, params, update=True):
    """Given a url and a tuple or dict of parameters, return
    a url that includes the parameters as a properly formatted
    query string.

    If update is True, change any existing values to new values
    given in params.
    """
    scheme, host, path, query, fragment = urlparse.urlsplit(url)

    # urlparse.parse_qsl gives back url-decoded byte strings. Leave these as
    # they are: they will be re-urlencoded below
    query_bits = [(k, v) for k, v in urlparse.parse_qsl(query)]
    if update:
        query_bits = dict(query_bits)
        query_bits.update(_make_query_tuples(params))
    else:
        query_bits.extend(_make_query_tuples(params))

    query = urllib.urlencode(query_bits)
    return urlparse.urlunsplit([scheme, host, path, query, fragment])


text:
# -*- coding: utf-8 -*-

import urllib2

import bunch
import ujson as json
from graph import GraphException
from url_operations import add_path, update_query_params

class FQL(object):

    """
    A maker of single and multiple FQL queries.

    Usage
    =====

    Single queries:

        >>> q = FQL('access_token')
        >>> result = q("SELECT post_id FROM stream WHERE source_id = ...")
        >>> result
        [Bunch(post_id='XXXYYYZZZ'), ...]

        >>> result[0]
        Bunch(post_id='XXXYYYZZZ')

        >>> result[0].post_id
        'XXXYYYZZZ'

    Multiple queries:

        >>> q = FQL('access_token')
        >>> result = q.multi(dict(query1="SELECT...", query2="SELECT..."))

        >>> result[0].name
        'query1'
        >>> result[0].fql_result_set
        [...]

        >>> result[1].name
        'query2'
        >>> result[1].fql_result_set
        [...]

    """

    ENDPOINT = 'https://api.facebook.com/method/'

    def __init__(self, access_token=None, err_handler=None):
        self.access_token = access_token
        self.err_handler = err_handler

    def __call__(self, query, **params):

        """
        Execute a single FQL query (using `fql.query`).

        Example:

            >>> q = FQL('access_token')
            >>> result = q("SELECT post_id FROM stream WHERE source_id = ...")
            >>> result
            [Bunch(post_id='XXXYYYZZZ'), ...]

            >>> result[0]
            Bunch(post_id='XXXYYYZZZ')

            >>> result[0].post_id
            'XXXYYYZZZ'

        """

        url = add_path(self.ENDPOINT, 'fql.query')
        params.update(query=query, access_token=self.access_token,
                      format='json')
        url = update_query_params(url, params)

        return self.fetch_json(url)

    def multi(self, queries, **params):

        """
        Execute multiple FQL queries (using `fql.multiquery`).

        Example:

            >>> q = FQL('access_token')
            >>> result = q.multi(dict(query1="SELECT...", query2="SELECT..."))

            >>> result[0].name
            'query1'
            >>> result[0].fql_result_set
            [...]

            >>> result[1].name
            'query2'
            >>> result[1].fql_result_set
            [...]

        """

        url = add_path(self.ENDPOINT, 'fql.multiquery')
        params.update(queries=json.dumps(queries),
                      access_token=self.access_token, format='json')
        url = update_query_params(url, params)

        return self.fetch_json(url)

    @classmethod
    def fetch_json(cls, url, data=None):
        response = json.loads(cls.fetch(url, data=data))
        if isinstance(response, dict):
            if response.get("error_msg"):
                code = response.get("error_code")
                msg = response.get("error_msg")
                args = response.get("request_args")
                raise GraphException(code, msg, args=args)
        return bunch.bunchify(response)

    @staticmethod
    def fetch(url, data=None):
        conn = urllib2.urlopen(url, data=data)
        try:
            return conn.read()
        finally:
            conn.close()