[med-svn] [ngila] 10/13: New upstream version 1.3
Andreas Tille
tille at debian.org
Wed Dec 13 15:45:44 UTC 2017
This is an automated email from the git hooks/post-receive script.
tille pushed a commit to branch master
in repository ngila.
commit b98b4024c61fbbba4ad9345dfc8bba8ef7549a94
Author: Andreas Tille <tille at debian.org>
Date: Wed Dec 13 16:41:18 2017 +0100
New upstream version 1.3
---
CMakeLists.txt | 218 +++++++++++
ChangeLog | 48 +++
Modules/FindGSL.cmake | 175 +++++++++
Modules/NSIS.template.in | 927 +++++++++++++++++++++++++++++++++++++++++++++++
copying.txt | 675 ++++++++++++++++++++++++++++++++++
debian/changelog | 5 -
debian/compat | 1 -
debian/control | 32 --
debian/copyright | 10 -
debian/rules | 10 -
debian/source/format | 1 -
debian/upstream/metadata | 12 -
debian/watch | 3 -
doc/CMakeLists.txt | 0
doc/mkman.pl | 41 +++
matrix/CMakeLists.txt | 0
matrix/dna | 5 +
ngilarc.txt | 4 +
readme.txt | 89 +++++
releng/build-rel-osx.sh | 55 +++
releng/build-rel-unix.sh | 45 +++
releng/build-rel-win.bat | 52 +++
src/CMakeLists.txt | 115 ++++++
src/align.cpp | 596 ++++++++++++++++++++++++++++++
src/align.h | 235 ++++++++++++
src/config.h.cmake | 10 +
src/lgmod.incl | 110 ++++++
src/matparser.cpp | 117 ++++++
src/matparser.h | 92 +++++
src/models.cpp | 442 ++++++++++++++++++++++
src/models.h | 116 ++++++
src/ngila-256.png | Bin 0 -> 34352 bytes
src/ngila-256.psd | Bin 0 -> 255054 bytes
src/ngila.cmds | 69 ++++
src/ngila.cpp | 437 ++++++++++++++++++++++
src/ngila.desktop | 13 +
src/ngila.h | 73 ++++
src/ngila.ico | Bin 0 -> 302430 bytes
src/ngila.png | Bin 0 -> 4425 bytes
src/ngila.rc | 1 +
src/ngila_app.h | 111 ++++++
src/seqdb.cpp | 147 ++++++++
src/seqdb.h | 190 ++++++++++
src/seqparser.h | 150 ++++++++
src/sort.h | 152 ++++++++
src/xm.h | 76 ++++
tests/CMakeLists.txt | 0
tests/test1.aln | 36 ++
tests/test1.fsa | 22 ++
tests/test1.sh | 10 +
tools/aamodel.R | 25 ++
tools/lg_LG.PAML.txt | 22 ++
tools/logzeta.R | 30 ++
53 files changed, 5731 insertions(+), 74 deletions(-)
diff --git a/CMakeLists.txt b/CMakeLists.txt
new file mode 100644
index 0000000..077cbe6
--- /dev/null
+++ b/CMakeLists.txt
@@ -0,0 +1,218 @@
+###########################################################
+# Ngila's Project Configuration
+#
+
+CMAKE_MINIMUM_REQUIRED(VERSION 2.6 FATAL_ERROR)
+SET(CMAKE_BUILD_TYPE_INIT "RelWithDebInfo")
+SET(CMAKE_VERBOSE_MAKEFILE OFF)
+SET(CMAKE_MODULE_PATH "${CMAKE_CURRENT_SOURCE_DIR}/Modules")
+
+IF(USE_STATIC_LIBS)
+ SET(Boost_USE_STATIC_LIBS ON)
+ SET(GSL_USE_STATIC_LIBS ON)
+ENDIF(USE_STATIC_LIBS)
+
+SET(USE_THREADS OFF CACHE BOOL
+ "Enable multi-threading support in Ngila (experimental).")
+MARK_AS_ADVANCED(USE_THREADS)
+
+PROJECT(Ngila)
+
+SET(PROJ_NAME "ngila")
+SET(PROJ_NAME_PRO "Ngila")
+SET(PROJ_EMAIL "reed at scit.us")
+SET(PROJ_VERSION_BRANCH release)
+SET(PROJ_VERSION_MAJOR 1)
+SET(PROJ_VERSION_MINOR 3)
+#SET(PROJ_VERSION_PATCH 0)
+SET(PROJ_VERSION_REV 0)
+
+IF(NOT PROJ_VERSION_REV)
+ # Test for Working copy
+ IF(EXISTS "${PROJECT_SOURCE_DIR}/.svn")
+ FIND_PACKAGE(Subversion)
+ IF(Subversion_FOUND)
+ Subversion_WC_INFO(${PROJECT_SOURCE_DIR} PROJ)
+ SET(PROJ_VERSION_REV "${PROJ_WC_LAST_CHANGED_REV}")
+ ENDIF(Subversion_FOUND)
+ ELSE(EXISTS "${PROJECT_SOURCE_DIR}/.svn")
+ # extract revision from source dir
+ STRING(REGEX REPLACE ".*-r([0-9]+)$" "\\1" PROJ_VERSION_REV ${PROJECT_SOURCE_DIR})
+ ENDIF(EXISTS "${PROJECT_SOURCE_DIR}/.svn")
+ENDIF(NOT PROJ_VERSION_REV)
+
+SET(PROJ_VERSION "${PROJ_VERSION_MAJOR}")
+IF(PROJ_VERSION_MINOR)
+ SET(PROJ_VERSION "${PROJ_VERSION}.${PROJ_VERSION_MINOR}")
+ENDIF(PROJ_VERSION_MINOR)
+IF(PROJ_VERSION_PATCH)
+ SET(PROJ_VERSION "${PROJ_VERSION}.${PROJ_VERSION_PATCH}")
+ENDIF(PROJ_VERSION_PATCH)
+IF(DEFINED PROJ_VERSION_BRANCH)
+ SET(PROJ_VERSION "${PROJ_VERSION}-${PROJ_VERSION_BRANCH}")
+ENDIF(DEFINED PROJ_VERSION_BRANCH)
+IF(PROJ_VERSION_REV)
+ IF(NOT PROJ_VERSION_BRANCH MATCHES release)
+ SET(PROJ_VERSION "${PROJ_VERSION}-r${PROJ_VERSION_REV}")
+ ENDIF(NOT PROJ_VERSION_BRANCH MATCHES release)
+ENDIF(PROJ_VERSION_REV)
+
+IF(APPLE)
+ IF(NOT DEFINED APPLE_BUNDLE)
+ SET(APPLE_BUNDLE ON)
+ ENDIF(NOT DEFINED APPLE_BUNDLE)
+ IF(NOT DEFINED APPLE_BUNDLE_APP)
+ SET(APPLE_BUNDLE_APP OFF)
+ ENDIF(NOT DEFINED APPLE_BUNDLE_APP)
+ELSE(APPLE)
+ SET(APPLE_BUNDLE OFF)
+ SET(APPLE_BUNDLE_APP OFF)
+ENDIF(APPLE)
+
+IF(APPLE_BUNDLE)
+ # Create 32-bit Universal Binary
+ # SET(CMAKE_OSX_ARCHITECTURES "ppc;i386")
+ # SET(CPACK_SYSTEM_NAME "${CMAKE_SYSTEM_NAME}-universal")
+ IF(CMAKE_INSTALL_PREFIX_INITIALIZED_TO_DEFAULT)
+ SET(CMAKE_INSTALL_PREFIX "/Applications")
+ ENDIF(CMAKE_INSTALL_PREFIX_INITIALIZED_TO_DEFAULT)
+ENDIF(APPLE_BUNDLE)
+
+IF(WIN32 AND NOT UNIX)
+ SET(PROJ_DIR)
+ELSEIF(APPLE_BUNDLE)
+ SET(PROJ_DIR)
+ELSE(WIN32 AND NOT UNIX)
+ SET(PROJ_DIR "/${PROJ_NAME}-${PROJ_VERSION}")
+ENDIF(WIN32 AND NOT UNIX)
+
+SET(CMAKE_DATA_DIR "share${PROJ_DIR}" CACHE STRING
+ "Install location for data (relative to prefix).")
+SET(CMAKE_DOC_DIR "doc${PROJ_DIR}" CACHE STRING
+ "Install location for documentation (relative to prefix).")
+SET(CMAKE_MAN_DIR "man" CACHE STRING
+ "Install location for man pages (relative to prefix).")
+MARK_AS_ADVANCED(CMAKE_DATA_DIR CMAKE_DOC_DIR CMAKE_MAN_DIR)
+
+SET(CPACK_SOURCE_IGNORE_FILES
+ "/CVS/" "/\\\\.svn/" "\\\\.swp$" "\\\\.#" "/#" ".*~$"
+ "/CMakeFiles/" "CMakeCache\\\\.txt"
+ "CPack.*Config\\\\.cmake" "cmake_install\\\\.cmake"
+ "install_manifest\\\\.txt$"
+ "_CPACK_PACKAGES" "_CPack_Packages"
+ "\\\\.vcproj" "\\\\.dir" "\\\\.ncb$" "\\\\.sln$" "\\\\.suo$"
+ "Makefile$" "\\\\.ilk" "\\\\.pdb"
+ "/releng/"
+ "${PROJ_NAME}[-]${PROJ_VERSION_MAJOR}" "${PROJ_NAME}$" "${PROJ_NAME}\\\\.exe"
+)
+
+SET(CPACK_PACKAGE_DESCRIPTION_SUMMARY "Pairwise alignments with log-affine gap costs")
+SET(CPACK_PACKAGE_VENDOR "Reed A. Cartwright")
+SET(CPACK_PACKAGE_EXECUTABLES "ngila" "Ngila")
+SET(CPACK_CREATE_DESKTOP_LINKS "ngila")
+
+SET(CPACK_PACKAGE_INSTALL_DIRECTORY ${PROJ_NAME_PRO})
+SET(CPACK_PACKAGE_NAME ${PROJ_NAME_PRO})
+SET(CPACK_PACKAGE_NAME_FILE ${PROJ_NAME})
+SET(CPACK_PACKAGE_VERSION ${PROJ_VERSION})
+SET(CPACK_PACKAGE_VERSION_MAJOR ${PROJ_VERSION_MAJOR})
+SET(CPACK_PACKAGE_VERSION_MINOR ${PROJ_VERSION_MINOR})
+SET(CPACK_PACKAGE_VERSION_PATCH ${PROJ_VERSION_PATCH})
+SET(CPACK_PACKAGE_DESCRIPTION_FILE "${CMAKE_CURRENT_SOURCE_DIR}/readme.txt")
+SET(CPACK_RESOURCE_FILE_LICENSE "${CMAKE_CURRENT_SOURCE_DIR}/copying.txt")
+
+SET(CPACK_SOURCE_PACKAGE_FILE_NAME "${CPACK_PACKAGE_NAME_FILE}-${CPACK_PACKAGE_VERSION}")
+
+IF(NOT DEFINED CPACK_SYSTEM_NAME)
+# make sure package is not Cygwin-unknown, for Cygwin just
+# cygwin is good for the system name
+ IF("${CMAKE_SYSTEM_NAME}" STREQUAL "CYGWIN")
+ SET(CPACK_SYSTEM_NAME Cygwin)
+ ELSE("${CMAKE_SYSTEM_NAME}" STREQUAL "CYGWIN")
+ SET(CPACK_SYSTEM_NAME ${CMAKE_SYSTEM_NAME}-${CMAKE_SYSTEM_PROCESSOR})
+ ENDIF("${CMAKE_SYSTEM_NAME}" STREQUAL "CYGWIN")
+ENDIF(NOT DEFINED CPACK_SYSTEM_NAME)
+IF(${CPACK_SYSTEM_NAME} MATCHES Windows)
+ IF(CMAKE_CL_64)
+ SET(CPACK_SYSTEM_NAME win64-${CMAKE_SYSTEM_PROCESSOR})
+ ELSE(CMAKE_CL_64)
+ SET(CPACK_SYSTEM_NAME win32-${CMAKE_SYSTEM_PROCESSOR})
+ ENDIF(CMAKE_CL_64)
+ENDIF(${CPACK_SYSTEM_NAME} MATCHES Windows)
+
+IF(NOT DEFINED CPACK_PACKAGE_FILE_NAME)
+# if the CPACK_PACKAGE_FILE_NAME is not defined by the cache
+# default to source package - system, on cygwin system is not
+# needed
+ IF(CYGWIN)
+ SET(CPACK_PACKAGE_FILE_NAME "${CPACK_SOURCE_PACKAGE_FILE_NAME}")
+ ELSE(CYGWIN)
+ SET(CPACK_PACKAGE_FILE_NAME
+ "${CPACK_SOURCE_PACKAGE_FILE_NAME}-${CPACK_SYSTEM_NAME}")
+ ENDIF(CYGWIN)
+ENDIF(NOT DEFINED CPACK_PACKAGE_FILE_NAME)
+
+IF(WIN32 AND NOT UNIX)
+ # There is a bug in NSI that does not handle full unix paths properly. Make
+ # sure there is at least one set of four (4) backlasshes.
+ #SET(CPACK_PACKAGE_ICON "${CMAKE_SOURCE_DIR}/Utilities/Release\\\\InstallIcon.bmp")
+ SET(CPACK_NSIS_INSTALLED_ICON_NAME "bin\\\\ngila.exe")
+ SET(CPACK_NSIS_HELP_LINK "http:\\\\\\\\scit.us\\\\projects\\\\ngila\\\\")
+ SET(CPACK_NSIS_URL_INFO_ABOUT "http:\\\\\\\\scit.us\\\\projects\\\\ngila\\\\")
+ SET(CPACK_NSIS_CONTACT "reed at scit.us")
+ SET(CPACK_NSIS_DISPLAY_NAME "${PROJ_NAME_PRO} ${PROJ_VERSION}")
+ SET(CPACK_NSIS_MODIFY_PATH ON)
+ SET(CPACK_CREATE_DESKTOP_LINKS)
+ SET(CPACK_NSIS_EXTRA_INSTALL_COMMANDS "
+ StrCmp \\\"$INSTALL_DESKTOP\\\" \\\"1\\\" 0 +2
+ CreateShortCut \\\"$DESKTOP\\\\Ngila.lnk\\\" \\\"$INSTDIR\\\\bin\\\\ngila.exe\\\" \\\"--desktop\\\"
+")
+ SET(CPACK_NSIS_EXTRA_UNINSTALL_COMMANDS "
+ StrCmp \\\"$INSTALL_DESKTOP\\\" \\\"1\\\" 0 +2
+ Delete \\\"$DESKTOP\\\\Ngila.lnk\\\"
+")
+
+ELSEIF(APPLE)
+ SET(CPACK_GENERATOR "PackageMaker")
+ SET(CPACK_SOURCE_GENERATOR "TGZ")
+ IF(APPLE_BUNDLE)
+ # make sure CMAKE_INSTALL_PREFIX ends in /
+ STRING(LENGTH "${CMAKE_INSTALL_PREFIX}" LEN)
+ MATH(EXPR LEN "${LEN} -1" )
+ STRING(SUBSTRING "${CMAKE_INSTALL_PREFIX}" ${LEN} 1 ENDCH)
+ IF(NOT "${ENDCH}" STREQUAL "/")
+ SET(CMAKE_INSTALL_PREFIX "${CMAKE_INSTALL_PREFIX}/")
+ ENDIF(NOT "${ENDCH}" STREQUAL "/")
+ SET(CPACK_PACKAGE_DEFAULT_LOCATION ${CMAKE_INSTALL_PREFIX})
+ SET(PROJ_BUNDLE_NAME "${PROJ_NAME} ${PROJ_VERSION}")
+ IF(APPLE_BUNDLE_APP)
+ SET(PROJ_BUNDLE_LOCATION "../..")
+ SET(CPACK_PACKAGING_INSTALL_PREFIX "/${PROJ_BUNDLE_NAME}.app/Contents")
+ SET(CMAKE_INSTALL_PREFIX
+ "${CMAKE_INSTALL_PREFIX}${CMAKE_BUNDLE_NAME}.app/Contents")
+ ELSE(APPLE_BUNDLE_APP)
+ SET(CPACK_PACKAGING_INSTALL_PREFIX "/${PROJ_BUNDLE_NAME}")
+ SET(CMAKE_INSTALL_PREFIX "${CMAKE_INSTALL_PREFIX}${CMAKE_BUNDLE_NAME}")
+ ENDIF(APPLE_BUNDLE_APP)
+ #CONFIGURE_FILE("postflight.sh.in" "postflight.sh")
+ #SET(CPACK_POSTFLIGHT_SCRIPT "postflight.sh")
+ ELSE(APPLE_BUNDLE)
+ ENDIF(APPLE_BUNDLE)
+ELSE(WIN32 AND NOT UNIX)
+ SET(CPACK_GENERATOR "TGZ")
+ SET(CPACK_SOURCE_GENERATOR "TGZ")
+ SET(CPACK_STRIP_FILES "bin/ngila")
+ENDIF(WIN32 AND NOT UNIX)
+
+INCLUDE(CPack)
+
+INSTALL(FILES copying.txt readme.txt ngilarc.txt ChangeLog DESTINATION ${CMAKE_DOC_DIR})
+#INSTALL(FILES lambda.pl outsplit.pl nexus2fasta.pl varrep.pl DESTINATION ${CMAKE_DATA_DIR})
+
+ADD_SUBDIRECTORY(src)
+ADD_SUBDIRECTORY(matrix)
+ADD_SUBDIRECTORY(doc)
+#ADD_SUBDIRECTORY(tests)
+
+#INCLUDE(InstallRequiredSystemLibraries)
+
diff --git a/ChangeLog b/ChangeLog
new file mode 100644
index 0000000..f285aca
--- /dev/null
+++ b/ChangeLog
@@ -0,0 +1,48 @@
+2010-12-09 1.3-RELEASE
+
+* Use CMake for compilation and installation
+* New scaling option enabled by default (identical sequences default to cost of 0)
+* Protein evolutionary models: aazeta and aageo
+* Fasta and Phylip format output support
+* Clustal and Phylip format input support
+* Report sequence identity measure
+* Matrix output formats for distance measures
+* Look for "ngilarc" file in the home directory.
+* New separator option
+* New const-align option
+* Replace arg-file option with ngilarc option.
+* Use custom zeta function if GSL not found.
+* Optimize size of travel table.
+
+* Ordering of --pairs-all fixed
+* bug fix for output of large alignments >10kb
+* minor bug fix for geo model
+
+
+2007-08-23 1.2.1-RELEASE
+
+* Fixed two bugs related to alignment printing
+
+2007-08-22 1.2-RELEASE
+
+* Changed license to GPL
+* Reorganized source code to a more modular framework.
+* Compiling from source now requires GSL and Boost Libraries
+* Added two evolutionary models: K2P+Zeta and K2P+Geometric
+* Command line options changed
+* Only FASTA supported for input
+* Tie breaking has changed
+* The Clustal output format has been altered
+* Man page updated and tied to cmds file
+* Added pairs setting for more control of sequential alignments
+* Several Bug fixes
+
+2007-02-23 1.1-RELEASE
+
+* Added options to specify free-end gap scores
+* Fixed bug with holisitc algorithm
+
+2006-10-26 1.0-RELEASE
+
+* Initial Release
+
diff --git a/Modules/FindGSL.cmake b/Modules/FindGSL.cmake
new file mode 100644
index 0000000..e505f52
--- /dev/null
+++ b/Modules/FindGSL.cmake
@@ -0,0 +1,175 @@
+# Try to find gnu scientific library GSL
+# See
+# http://www.gnu.org/software/gsl/ and
+# http://gnuwin32.sourceforge.net/packages/gsl.htm
+#
+# Once run this will define:
+#
+# GSL_FOUND = system has GSL lib
+#
+# GSL_LIBRARIES = full path to the libraries
+# on Unix/Linux with additional linker flags from "gsl-config --libs"
+#
+# CMAKE_GSL_C_FLAGS = Unix compiler flags for GSL, essentially "`gsl-config --cflags`"
+#
+# GSL_INCLUDE_DIR = where to find headers
+#
+# GSL_LINK_DIRECTORIES = link directories, useful for rpath on Unix
+# GSL_EXE_LINKER_FLAGS = rpath on Unix
+#
+# Felix Woelk 07/2004
+# Jan Woetzel
+#
+# www.mip.informatik.uni-kiel.de
+# --------------------------------
+
+IF(WIN32)
+ # JW tested with gsl-1.8, Windows XP, MSVS 7.1
+ SET(GSL_POSSIBLE_ROOT_DIRS
+ ${GSL_ROOT_DIR}
+ $ENV{GSL_ROOT_DIR}
+ ${GSL_DIR}
+ ${GSL_HOME}
+ $ENV{GSL_DIR}
+ $ENV{GSL_HOME}
+ $ENV{EXTRA}
+ "C:/Programs/gsl/"
+ "C:/Program Files (x86)/gsl/"
+ "C:/Program Files/gsl/"
+ "C:/Programs/"
+ "C:/Program Files/"
+ )
+ IF(GSL_USE_STATIC_LIBS)
+ SET(GSLSTATIC lib)
+ ENDIF(GSL_USE_STATIC_LIBS)
+
+
+ FIND_PATH(GSL_INCLUDE_DIR
+ NAMES gsl/gsl_cdf.h gsl/gsl_randist.h
+ HINTS ${GSL_POSSIBLE_ROOT_DIRS}
+ PATH_SUFFIXES include
+ DOC "GSL header include dir"
+ )
+
+ FIND_LIBRARY(GSL_GSL_LIBRARY
+ NAMES "${GSLSTATIC}gsl"
+ HINTS ${GSL_POSSIBLE_ROOT_DIRS}
+ PATH_SUFFIXES lib
+ DOC "GSL library dir" )
+
+ FIND_LIBRARY(GSL_GSLCBLAS_LIBRARY
+ NAMES "${GSLSTATIC}gslcblas" "${GSLSTATIC}cblas"
+ HINTS ${GSL_POSSIBLE_ROOT_DIRS}
+ PATH_SUFFIXES lib
+ DOC "GSL cblas library dir" )
+
+ SET(GSL_LIBRARIES ${GSL_GSL_LIBRARY})
+
+ MARK_AS_ADVANCED(
+ GSL_GSLCBLAS_LIBRARY
+ GSL_GSL_LIBRARY
+ GSL_INCLUDE_DIR
+ )
+
+
+ELSE(WIN32)
+
+ IF(UNIX)
+ SET(GSL_CONFIG_PREFER_PATH
+ "$ENV{GSL_DIR}/bin"
+ "$ENV{GSL_DIR}"
+ "$ENV{GSL_HOME}/bin"
+ "$ENV{GSL_HOME}"
+ CACHE STRING "preferred path to GSL (gsl-config)")
+ FIND_PROGRAM(GSL_CONFIG gsl-config
+ HINTS ${GSL_CONFIG_PREFER_PATH}
+ /usr/bin/
+ DOC "Path to gsl-config binary"
+ )
+ # MESSAGE("DBG GSL_CONFIG ${GSL_CONFIG}")
+
+ IF (GSL_CONFIG)
+ # set CFLAGS to be fed into C_FLAGS by the user:
+ EXEC_PROGRAM(${GSL_CONFIG}
+ ARGS --cflags
+ OUTPUT_VARIABLE GSL_CONFIG_CFLAGS )
+ SET(GSL_C_FLAGS ${GSL_CONFIG_CFLAGS}
+ CACHE STRING "Flags used to compile against GSL")
+
+ IF (NOT GSL_PREFIX)
+ EXEC_PROGRAM(${GSL_CONFIG}
+ ARGS --prefix
+ OUTPUT_VARIABLE GSL_PREFIX)
+ SET(GSL_PREFIX ${GSL_PREFIX} CACHE STRING "Location of GSL")
+ MESSAGE(STATUS "Using GSL from ${GSL_PREFIX}")
+ ENDIF(NOT GSL_PREFIX)
+
+ # set INCLUDE_DIRS to prefix+include
+ SET(GSL_INCLUDE_DIR ${GSL_PREFIX}/include
+ CACHE STRING "Location of GSL headers")
+
+
+ # set link libraries and link flags
+ EXEC_PROGRAM(${GSL_CONFIG}
+ ARGS --libs
+ OUTPUT_VARIABLE GSL_CONFIG_LIBS )
+ STRING(REGEX MATCHALL "[-][l][^ ;]+"
+ GSL_LIBRARIES_WITH_PREFIX
+ "${GSL_CONFIG_LIBS}" )
+ IF(GSL_USE_STATIC_LIBS)
+ IF(UNIX)
+ STRING(REGEX REPLACE "[-][l]([^ ;]+)" "lib\\1.a" GSL_LIBRARIES_WITH_PREFIX "${GSL_LIBRARIES_WITH_PREFIX}")
+ ENDIF(UNIX)
+ ENDIF(GSL_USE_STATIC_LIBS)
+ SET(GSL_LIBRARIES ${GSL_LIBRARIES_WITH_PREFIX}
+ CACHE STRING "GSL libraries to link against")
+
+ # extract link dirs for rpath
+ # split off the link dirs (for rpath)
+ # use regular expression to match wildcard equivalent "-L*<endchar>"
+ # with <endchar> is a space or a semicolon
+ STRING(REGEX MATCHALL "[-][L]([^ ;])+"
+ GSL_LINK_DIRECTORIES_WITH_PREFIX
+ "${GSL_CONFIG_LIBS}" )
+ # MESSAGE("DBG GSL_LINK_DIRECTORIES_WITH_PREFIX=${GSL_LINK_DIRECTORIES_WITH_PREFIX}")
+
+ # remove prefix -L because we need the pure directory for LINK_DIRECTORIES
+
+ IF (GSL_LINK_DIRECTORIES_WITH_PREFIX)
+ STRING(REGEX REPLACE "[-][L]" "" GSL_LINK_DIRECTORIES ${GSL_LINK_DIRECTORIES_WITH_PREFIX} )
+ ENDIF (GSL_LINK_DIRECTORIES_WITH_PREFIX)
+ SET(GSL_EXE_LINKER_FLAGS "-Wl,-rpath,${GSL_LINK_DIRECTORIES}"
+ CACHE STRING "Flags used when linking against GSL")
+ # MESSAGE("DBG GSL_LINK_DIRECTORIES=${GSL_LINK_DIRECTORIES}")
+ # MESSAGE("DBG GSL_EXE_LINKER_FLAGS=${GSL_EXE_LINKER_FLAGS}")
+
+ # ADD_DEFINITIONS("-DHAVE_GSL")
+ # SET(GSL_DEFINITIONS "-DHAVE_GSL")
+ MARK_AS_ADVANCED(
+ GSL_C_FLAGS
+ GSL_INCLUDE_DIR
+ GSL_LIBRARIES
+ GSL_LINK_DIRECTORIES
+ GSL_DEFINITIONS
+ GSL_EXE_LINKER_FLAGS
+ GSL_CONFIG
+ GSL_PREFIX
+ GSL_CONFIG_PREFER_PATH
+ )
+
+ ELSE(GSL_CONFIG)
+ MESSAGE("FindGSL.cmake: gsl-config not found. Please set it manually. GSL_CONFIG=${GSL_CONFIG}")
+ ENDIF(GSL_CONFIG)
+
+ ENDIF(UNIX)
+ENDIF(WIN32)
+
+
+IF(GSL_LIBRARIES)
+ IF(GSL_INCLUDE_DIR OR GSL_C_FLAGS)
+
+ SET(GSL_FOUND 1)
+
+ ENDIF(GSL_INCLUDE_DIR OR GSL_C_FLAGS)
+ENDIF(GSL_LIBRARIES)
+
diff --git a/Modules/NSIS.template.in b/Modules/NSIS.template.in
new file mode 100644
index 0000000..56e59f6
--- /dev/null
+++ b/Modules/NSIS.template.in
@@ -0,0 +1,927 @@
+; CPack install script designed for a nmake build
+
+;--------------------------------
+; You must define these values
+
+ !define VERSION "@CPACK_PACKAGE_VERSION@"
+ !define PATCH "@CPACK_PACKAGE_VERSION_PATCH@"
+ !define INST_DIR "@CPACK_TEMPORARY_DIRECTORY@"
+
+;--------------------------------
+;Variables
+
+ Var MUI_TEMP
+ Var STARTMENU_FOLDER
+ Var SV_ALLUSERS
+ Var START_MENU
+ Var DO_NOT_ADD_TO_PATH
+ Var ADD_TO_PATH_ALL_USERS
+ Var ADD_TO_PATH_CURRENT_USER
+ Var INSTALL_DESKTOP
+
+;-------------------------------
+;Vista Orange Theme
+;Sets the theme path
+ !define OMUI_THEME_PATH "${NSISDIR}\Contrib\MUI Orange Vista Theme\Clean"
+
+; MUI Settings / Icons
+; In the moment of writing this, NSIS don't support well Vista icons with PNG compression.
+; We provide both, compressed and uncompressed (-nopng) icons.
+ !define MUI_ICON "${OMUI_THEME_PATH}\installer-nopng.ico"
+ !define MUI_UNICON "${OMUI_THEME_PATH}\uninstaller-nopng.ico"
+
+; MUI Settings / Header
+ !define MUI_HEADERIMAGE
+ !define MUI_HEADERIMAGE_RIGHT
+ !define MUI_HEADERIMAGE_BITMAP "${OMUI_THEME_PATH}\header-r.bmp"
+ !define MUI_HEADERIMAGE_UNBITMAP "${OMUI_THEME_PATH}\header-r-un.bmp"
+
+; MUI Settings / Wizard
+ !define MUI_WELCOMEFINISHPAGE_BITMAP "${OMUI_THEME_PATH}\wizard.bmp"
+ !define MUI_UNWELCOMEFINISHPAGE_BITMAP "${OMUI_THEME_PATH}\wizard-un.bmp"
+
+;--------------------------------
+;Include Modern UI
+
+ !include "MUI.nsh"
+
+ ;Default installation folder
+ InstallDir "$PROGRAMFILES\@CPACK_PACKAGE_INSTALL_DIRECTORY@"
+
+;--------------------------------
+;General
+
+ ;Name and file
+ Name "@CPACK_PACKAGE_INSTALL_DIRECTORY@"
+ OutFile "@CPACK_TOPLEVEL_DIRECTORY@/@CPACK_OUTPUT_FILE_NAME@"
+
+ ;Set compression
+ SetCompressor @CPACK_NSIS_COMPRESSOR@
+
+ at CPACK_NSIS_DEFINES@
+
+ !include Sections.nsh
+
+;--- Component support macros: ---
+; The code for the add/remove functionality is from:
+; http://nsis.sourceforge.net/Add/Remove_Functionality
+; It has been modified slightly and extended to provide
+; inter-component dependencies.
+Var AR_SecFlags
+Var AR_RegFlags
+ at CPACK_NSIS_SECTION_SELECTED_VARS@
+
+; Loads the "selected" flag for the section named SecName into the
+; variable VarName.
+!macro LoadSectionSelectedIntoVar SecName VarName
+ SectionGetFlags ${${SecName}} $${VarName}
+ IntOp $${VarName} $${VarName} & ${SF_SELECTED} ;Turn off all other bits
+!macroend
+
+; Loads the value of a variable... can we get around this?
+!macro LoadVar VarName
+ IntOp $R0 0 + $${VarName}
+!macroend
+
+; Sets the value of a variable
+!macro StoreVar VarName IntValue
+ IntOp $${VarName} 0 + ${IntValue}
+!macroend
+
+!macro InitSection SecName
+ ; This macro reads component installed flag from the registry and
+ ;changes checked state of the section on the components page.
+ ;Input: section index constant name specified in Section command.
+
+ ClearErrors
+ ;Reading component status from registry
+ ReadRegDWORD $AR_RegFlags HKLM "Software\Microsoft\Windows\CurrentVersion\Uninstall\@CPACK_PACKAGE_INSTALL_DIRECTORY@\Components\${SecName}" "Installed"
+ IfErrors "default_${SecName}"
+ ;Status will stay default if registry value not found
+ ;(component was never installed)
+ IntOp $AR_RegFlags $AR_RegFlags & ${SF_SELECTED} ;Turn off all other bits
+ SectionGetFlags ${${SecName}} $AR_SecFlags ;Reading default section flags
+ IntOp $AR_SecFlags $AR_SecFlags & 0xFFFE ;Turn lowest (enabled) bit off
+ IntOp $AR_SecFlags $AR_RegFlags | $AR_SecFlags ;Change lowest bit
+
+ ; Note whether this component was installed before
+ !insertmacro StoreVar ${SecName}_was_installed $AR_RegFlags
+ IntOp $R0 $AR_RegFlags & $AR_RegFlags
+
+ ;Writing modified flags
+ SectionSetFlags ${${SecName}} $AR_SecFlags
+
+ "default_${SecName}:"
+ !insertmacro LoadSectionSelectedIntoVar ${SecName} ${SecName}_selected
+!macroend
+
+!macro FinishSection SecName
+ ; This macro reads section flag set by user and removes the section
+ ;if it is not selected.
+ ;Then it writes component installed flag to registry
+ ;Input: section index constant name specified in Section command.
+
+ SectionGetFlags ${${SecName}} $AR_SecFlags ;Reading section flags
+ ;Checking lowest bit:
+ IntOp $AR_SecFlags $AR_SecFlags & ${SF_SELECTED}
+ IntCmp $AR_SecFlags 1 "leave_${SecName}"
+ ;Section is not selected:
+ ;Calling Section uninstall macro and writing zero installed flag
+ !insertmacro "Remove_${${SecName}}"
+ WriteRegDWORD HKLM "Software\Microsoft\Windows\CurrentVersion\Uninstall\@CPACK_PACKAGE_INSTALL_DIRECTORY@\Components\${SecName}" \
+ "Installed" 0
+ Goto "exit_${SecName}"
+
+ "leave_${SecName}:"
+ ;Section is selected:
+ WriteRegDWORD HKLM "Software\Microsoft\Windows\CurrentVersion\Uninstall\@CPACK_PACKAGE_INSTALL_DIRECTORY@\Components\${SecName}" \
+ "Installed" 1
+
+ "exit_${SecName}:"
+!macroend
+
+!macro RemoveSection SecName
+ ; This macro is used to call section's Remove_... macro
+ ;from the uninstaller.
+ ;Input: section index constant name specified in Section command.
+
+ !insertmacro "Remove_${${SecName}}"
+!macroend
+
+; Determine whether the selection of SecName changed
+!macro MaybeSelectionChanged SecName
+ !insertmacro LoadVar ${SecName}_selected
+ SectionGetFlags ${${SecName}} $R1
+ IntOp $R1 $R1 & ${SF_SELECTED} ;Turn off all other bits
+
+ ; See if the status has changed:
+ IntCmp $R0 $R1 "${SecName}_unchanged"
+ !insertmacro LoadSectionSelectedIntoVar ${SecName} ${SecName}_selected
+
+ IntCmp $R1 ${SF_SELECTED} "${SecName}_was_selected"
+ !insertmacro "Deselect_required_by_${SecName}"
+ goto "${SecName}_unchanged"
+
+ "${SecName}_was_selected:"
+ !insertmacro "Select_${SecName}_depends"
+
+ "${SecName}_unchanged:"
+!macroend
+;--- End of Add/Remove macros ---
+
+;--------------------------------
+;Interface Settings
+
+; !define MUI_HEADERIMAGE
+ !define MUI_ABORTWARNING
+
+;--------------------------------
+; path functions
+
+!verbose 3
+!include "WinMessages.NSH"
+!verbose 4
+
+;----------------------------------------
+; based upon a script of "Written by KiCHiK 2003-01-18 05:57:02"
+;----------------------------------------
+!verbose 3
+!include "WinMessages.NSH"
+!verbose 4
+;====================================================
+; get_NT_environment
+; Returns: the selected environment
+; Output : head of the stack
+;====================================================
+!macro select_NT_profile UN
+Function ${UN}select_NT_profile
+ StrCmp $ADD_TO_PATH_ALL_USERS "1" 0 environment_single
+ DetailPrint "Selected environment for all users"
+ Push "all"
+ Return
+ environment_single:
+ DetailPrint "Selected environment for current user only."
+ Push "current"
+ Return
+FunctionEnd
+!macroend
+!insertmacro select_NT_profile ""
+!insertmacro select_NT_profile "un."
+;----------------------------------------------------
+!define NT_current_env 'HKCU "Environment"'
+!define NT_all_env 'HKLM "SYSTEM\CurrentControlSet\Control\Session Manager\Environment"'
+
+!ifndef WriteEnvStr_RegKey
+ !ifdef ALL_USERS
+ !define WriteEnvStr_RegKey \
+ 'HKLM "SYSTEM\CurrentControlSet\Control\Session Manager\Environment"'
+ !else
+ !define WriteEnvStr_RegKey 'HKCU "Environment"'
+ !endif
+!endif
+
+; AddToPath - Adds the given dir to the search path.
+; Input - head of the stack
+; Note - Win9x systems requires reboot
+
+Function AddToPath
+ Exch $0
+ Push $1
+ Push $2
+ Push $3
+
+ # don't add if the path doesn't exist
+ IfFileExists "$0\*.*" "" AddToPath_done
+
+ ReadEnvStr $1 PATH
+ Push "$1;"
+ Push "$0;"
+ Call StrStr
+ Pop $2
+ StrCmp $2 "" "" AddToPath_done
+ Push "$1;"
+ Push "$0\;"
+ Call StrStr
+ Pop $2
+ StrCmp $2 "" "" AddToPath_done
+ GetFullPathName /SHORT $3 $0
+ Push "$1;"
+ Push "$3;"
+ Call StrStr
+ Pop $2
+ StrCmp $2 "" "" AddToPath_done
+ Push "$1;"
+ Push "$3\;"
+ Call StrStr
+ Pop $2
+ StrCmp $2 "" "" AddToPath_done
+
+ Call IsNT
+ Pop $1
+ StrCmp $1 1 AddToPath_NT
+ ; Not on NT
+ StrCpy $1 $WINDIR 2
+ FileOpen $1 "$1\autoexec.bat" a
+ FileSeek $1 -1 END
+ FileReadByte $1 $2
+ IntCmp $2 26 0 +2 +2 # DOS EOF
+ FileSeek $1 -1 END # write over EOF
+ FileWrite $1 "$\r$\nSET PATH=%PATH%;$3$\r$\n"
+ FileClose $1
+ SetRebootFlag true
+ Goto AddToPath_done
+
+ AddToPath_NT:
+ ReadRegStr $1 ${WriteEnvStr_RegKey} "PATH"
+ StrCmp $1 "" AddToPath_NTdoIt
+ Push $1
+ Call Trim
+ Pop $1
+ StrCpy $0 "$1;$0"
+ AddToPath_NTdoIt:
+ WriteRegExpandStr ${WriteEnvStr_RegKey} "PATH" $0
+ SendMessage ${HWND_BROADCAST} ${WM_WININICHANGE} 0 "STR:Environment" /TIMEOUT=5000
+
+ AddToPath_done:
+ Pop $3
+ Pop $2
+ Pop $1
+ Pop $0
+FunctionEnd
+
+
+; RemoveFromPath - Remove a given dir from the path
+; Input: head of the stack
+
+Function un.RemoveFromPath
+ Exch $0
+ Push $1
+ Push $2
+ Push $3
+ Push $4
+ Push $5
+ Push $6
+
+ IntFmt $6 "%c" 26 # DOS EOF
+
+ Call un.IsNT
+ Pop $1
+ StrCmp $1 1 unRemoveFromPath_NT
+ ; Not on NT
+ StrCpy $1 $WINDIR 2
+ FileOpen $1 "$1\autoexec.bat" r
+ GetTempFileName $4
+ FileOpen $2 $4 w
+ GetFullPathName /SHORT $0 $0
+ StrCpy $0 "SET PATH=%PATH%;$0"
+ Goto unRemoveFromPath_dosLoop
+
+ unRemoveFromPath_dosLoop:
+ FileRead $1 $3
+ StrCpy $5 $3 1 -1 # read last char
+ StrCmp $5 $6 0 +2 # if DOS EOF
+ StrCpy $3 $3 -1 # remove DOS EOF so we can compare
+ StrCmp $3 "$0$\r$\n" unRemoveFromPath_dosLoopRemoveLine
+ StrCmp $3 "$0$\n" unRemoveFromPath_dosLoopRemoveLine
+ StrCmp $3 "$0" unRemoveFromPath_dosLoopRemoveLine
+ StrCmp $3 "" unRemoveFromPath_dosLoopEnd
+ FileWrite $2 $3
+ Goto unRemoveFromPath_dosLoop
+ unRemoveFromPath_dosLoopRemoveLine:
+ SetRebootFlag true
+ Goto unRemoveFromPath_dosLoop
+
+ unRemoveFromPath_dosLoopEnd:
+ FileClose $2
+ FileClose $1
+ StrCpy $1 $WINDIR 2
+ Delete "$1\autoexec.bat"
+ CopyFiles /SILENT $4 "$1\autoexec.bat"
+ Delete $4
+ Goto unRemoveFromPath_done
+
+ unRemoveFromPath_NT:
+ ReadRegStr $1 ${WriteEnvStr_RegKey} "PATH"
+ StrCpy $5 $1 1 -1 # copy last char
+ StrCmp $5 ";" +2 # if last char != ;
+ StrCpy $1 "$1;" # append ;
+ Push $1
+ Push "$0;"
+ Call un.StrStr ; Find `$0;` in $1
+ Pop $2 ; pos of our dir
+ StrCmp $2 "" unRemoveFromPath_done
+ ; else, it is in path
+ # $0 - path to add
+ # $1 - path var
+ StrLen $3 "$0;"
+ StrLen $4 $2
+ StrCpy $5 $1 -$4 # $5 is now the part before the path to remove
+ StrCpy $6 $2 "" $3 # $6 is now the part after the path to remove
+ StrCpy $3 $5$6
+
+ StrCpy $5 $3 1 -1 # copy last char
+ StrCmp $5 ";" 0 +2 # if last char == ;
+ StrCpy $3 $3 -1 # remove last char
+
+ WriteRegExpandStr ${WriteEnvStr_RegKey} "PATH" $3
+ SendMessage ${HWND_BROADCAST} ${WM_WININICHANGE} 0 "STR:Environment" /TIMEOUT=5000
+
+ unRemoveFromPath_done:
+ Pop $6
+ Pop $5
+ Pop $4
+ Pop $3
+ Pop $2
+ Pop $1
+ Pop $0
+FunctionEnd
+
+;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+; Uninstall sutff
+;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+
+###########################################
+# Utility Functions #
+###########################################
+
+;====================================================
+; IsNT - Returns 1 if the current system is NT, 0
+; otherwise.
+; Output: head of the stack
+;====================================================
+; IsNT
+; no input
+; output, top of the stack = 1 if NT or 0 if not
+;
+; Usage:
+; Call IsNT
+; Pop $R0
+; ($R0 at this point is 1 or 0)
+
+!macro IsNT un
+Function ${un}IsNT
+ Push $0
+ ReadRegStr $0 HKLM "SOFTWARE\Microsoft\Windows NT\CurrentVersion" CurrentVersion
+ StrCmp $0 "" 0 IsNT_yes
+ ; we are not NT.
+ Pop $0
+ Push 0
+ Return
+
+ IsNT_yes:
+ ; NT!!!
+ Pop $0
+ Push 1
+FunctionEnd
+!macroend
+!insertmacro IsNT ""
+!insertmacro IsNT "un."
+
+; StrStr
+; input, top of stack = string to search for
+; top of stack-1 = string to search in
+; output, top of stack (replaces with the portion of the string remaining)
+; modifies no other variables.
+;
+; Usage:
+; Push "this is a long ass string"
+; Push "ass"
+; Call StrStr
+; Pop $R0
+; ($R0 at this point is "ass string")
+
+!macro StrStr un
+Function ${un}StrStr
+Exch $R1 ; st=haystack,old$R1, $R1=needle
+ Exch ; st=old$R1,haystack
+ Exch $R2 ; st=old$R1,old$R2, $R2=haystack
+ Push $R3
+ Push $R4
+ Push $R5
+ StrLen $R3 $R1
+ StrCpy $R4 0
+ ; $R1=needle
+ ; $R2=haystack
+ ; $R3=len(needle)
+ ; $R4=cnt
+ ; $R5=tmp
+ loop:
+ StrCpy $R5 $R2 $R3 $R4
+ StrCmp $R5 $R1 done
+ StrCmp $R5 "" done
+ IntOp $R4 $R4 + 1
+ Goto loop
+done:
+ StrCpy $R1 $R2 "" $R4
+ Pop $R5
+ Pop $R4
+ Pop $R3
+ Pop $R2
+ Exch $R1
+FunctionEnd
+!macroend
+!insertmacro StrStr ""
+!insertmacro StrStr "un."
+
+Function Trim ; Added by Pelaca
+ Exch $R1
+ Push $R2
+Loop:
+ StrCpy $R2 "$R1" 1 -1
+ StrCmp "$R2" " " RTrim
+ StrCmp "$R2" "$\n" RTrim
+ StrCmp "$R2" "$\r" RTrim
+ StrCmp "$R2" ";" RTrim
+ GoTo Done
+RTrim:
+ StrCpy $R1 "$R1" -1
+ Goto Loop
+Done:
+ Pop $R2
+ Exch $R1
+FunctionEnd
+
+Function ConditionalAddToRegisty
+ Pop $0
+ Pop $1
+ StrCmp "$0" "" ConditionalAddToRegisty_EmptyString
+ WriteRegStr SHCTX "Software\Microsoft\Windows\CurrentVersion\Uninstall\@CPACK_PACKAGE_INSTALL_DIRECTORY@" \
+ "$1" "$0"
+ ;MessageBox MB_OK "Set Registry: '$1' to '$0'"
+ DetailPrint "Set install registry entry: '$1' to '$0'"
+ ConditionalAddToRegisty_EmptyString:
+FunctionEnd
+
+;--------------------------------
+
+!ifdef CPACK_USES_DOWNLOAD
+Function DownloadFile
+ IfFileExists $INSTDIR\* +2
+ CreateDirectory $INSTDIR
+ Pop $0
+
+ ; Skip if already downloaded
+ IfFileExists $INSTDIR\$0 0 +2
+ Return
+
+ StrCpy $1 "@CPACK_DOWNLOAD_SITE@"
+
+ try_again:
+ NSISdl::download "$1/$0" "$INSTDIR\$0"
+
+ Pop $1
+ StrCmp $1 "success" success
+ StrCmp $1 "Cancelled" cancel
+ MessageBox MB_OK "Download failed: $1"
+ cancel:
+ Return
+ success:
+FunctionEnd
+!endif
+
+;--------------------------------
+; Installation types
+ at CPACK_NSIS_INSTALLATION_TYPES@
+
+;--------------------------------
+; Component sections
+ at CPACK_NSIS_COMPONENT_SECTIONS@
+
+;--------------------------------
+; Define some macro setting for the gui
+ at CPACK_NSIS_INSTALLER_MUI_ICON_CODE@
+ at CPACK_NSIS_INSTALLER_ICON_CODE@
+ at CPACK_NSIS_INSTALLER_MUI_COMPONENTS_DESC@
+
+;--------------------------------
+;Pages
+ !insertmacro MUI_PAGE_WELCOME
+
+ !insertmacro MUI_PAGE_LICENSE "@CPACK_RESOURCE_FILE_LICENSE@"
+ Page custom InstallOptionsPage
+ !insertmacro MUI_PAGE_DIRECTORY
+
+ ;Start Menu Folder Page Configuration
+ !define MUI_STARTMENUPAGE_REGISTRY_ROOT "SHCTX"
+ !define MUI_STARTMENUPAGE_REGISTRY_KEY "Software\@CPACK_PACKAGE_VENDOR@\@CPACK_PACKAGE_INSTALL_REGISTRY_KEY@"
+ !define MUI_STARTMENUPAGE_REGISTRY_VALUENAME "Start Menu Folder"
+ !insertmacro MUI_PAGE_STARTMENU Application $STARTMENU_FOLDER
+
+ @CPACK_NSIS_PAGE_COMPONENTS@
+
+ !insertmacro MUI_PAGE_INSTFILES
+ !insertmacro MUI_PAGE_FINISH
+
+ !insertmacro MUI_UNPAGE_CONFIRM
+ !insertmacro MUI_UNPAGE_INSTFILES
+
+;--------------------------------
+;Languages
+
+ !insertmacro MUI_LANGUAGE "English" ;first language is the default language
+ !insertmacro MUI_LANGUAGE "Albanian"
+ !insertmacro MUI_LANGUAGE "Arabic"
+ !insertmacro MUI_LANGUAGE "Basque"
+ !insertmacro MUI_LANGUAGE "Belarusian"
+ !insertmacro MUI_LANGUAGE "Bosnian"
+ !insertmacro MUI_LANGUAGE "Breton"
+ !insertmacro MUI_LANGUAGE "Bulgarian"
+ !insertmacro MUI_LANGUAGE "Croatian"
+ !insertmacro MUI_LANGUAGE "Czech"
+ !insertmacro MUI_LANGUAGE "Danish"
+ !insertmacro MUI_LANGUAGE "Dutch"
+ !insertmacro MUI_LANGUAGE "Estonian"
+ !insertmacro MUI_LANGUAGE "Farsi"
+ !insertmacro MUI_LANGUAGE "Finnish"
+ !insertmacro MUI_LANGUAGE "French"
+ !insertmacro MUI_LANGUAGE "German"
+ !insertmacro MUI_LANGUAGE "Greek"
+ !insertmacro MUI_LANGUAGE "Hebrew"
+ !insertmacro MUI_LANGUAGE "Hungarian"
+ !insertmacro MUI_LANGUAGE "Icelandic"
+ !insertmacro MUI_LANGUAGE "Indonesian"
+ !insertmacro MUI_LANGUAGE "Irish"
+ !insertmacro MUI_LANGUAGE "Italian"
+ !insertmacro MUI_LANGUAGE "Japanese"
+ !insertmacro MUI_LANGUAGE "Korean"
+ !insertmacro MUI_LANGUAGE "Kurdish"
+ !insertmacro MUI_LANGUAGE "Latvian"
+ !insertmacro MUI_LANGUAGE "Lithuanian"
+ !insertmacro MUI_LANGUAGE "Luxembourgish"
+ !insertmacro MUI_LANGUAGE "Macedonian"
+ !insertmacro MUI_LANGUAGE "Malay"
+ !insertmacro MUI_LANGUAGE "Mongolian"
+ !insertmacro MUI_LANGUAGE "Norwegian"
+ !insertmacro MUI_LANGUAGE "Polish"
+ !insertmacro MUI_LANGUAGE "Portuguese"
+ !insertmacro MUI_LANGUAGE "PortugueseBR"
+ !insertmacro MUI_LANGUAGE "Romanian"
+ !insertmacro MUI_LANGUAGE "Russian"
+ !insertmacro MUI_LANGUAGE "Serbian"
+ !insertmacro MUI_LANGUAGE "SerbianLatin"
+ !insertmacro MUI_LANGUAGE "SimpChinese"
+ !insertmacro MUI_LANGUAGE "Slovak"
+ !insertmacro MUI_LANGUAGE "Slovenian"
+ !insertmacro MUI_LANGUAGE "Spanish"
+ !insertmacro MUI_LANGUAGE "Swedish"
+ !insertmacro MUI_LANGUAGE "Thai"
+ !insertmacro MUI_LANGUAGE "TradChinese"
+ !insertmacro MUI_LANGUAGE "Turkish"
+ !insertmacro MUI_LANGUAGE "Ukrainian"
+ !insertmacro MUI_LANGUAGE "Welsh"
+
+
+;--------------------------------
+;Reserve Files
+
+ ;These files should be inserted before other files in the data block
+ ;Keep these lines before any File command
+ ;Only for solid compression (by default, solid compression is enabled for BZIP2 and LZMA)
+
+ ReserveFile "NSIS.InstallOptions.ini"
+ !insertmacro MUI_RESERVEFILE_INSTALLOPTIONS
+
+;--------------------------------
+;Installer Sections
+
+Section "-Core installation"
+ ;Use the entire tree produced by the INSTALL target. Keep the
+ ;list of directories here in sync with the RMDir commands below.
+ SetOutPath "$INSTDIR"
+ @CPACK_NSIS_FULL_INSTALL@
+
+ ;Store installation folder
+ WriteRegStr SHCTX "Software\@CPACK_PACKAGE_VENDOR@\@CPACK_PACKAGE_INSTALL_REGISTRY_KEY@" "" $INSTDIR
+
+ ;Create uninstaller
+ WriteUninstaller "$INSTDIR\Uninstall.exe"
+ Push "DisplayName"
+ Push "@CPACK_NSIS_DISPLAY_NAME@"
+ Call ConditionalAddToRegisty
+ Push "DisplayVersion"
+ Push "@CPACK_PACKAGE_VERSION@"
+ Call ConditionalAddToRegisty
+ Push "Publisher"
+ Push "@CPACK_PACKAGE_VENDOR@"
+ Call ConditionalAddToRegisty
+ Push "UninstallString"
+ Push "$INSTDIR\Uninstall.exe"
+ Call ConditionalAddToRegisty
+ Push "NoRepair"
+ Push "1"
+ Call ConditionalAddToRegisty
+
+ !ifdef CPACK_NSIS_ADD_REMOVE
+ ;Create add/remove functionality
+ Push "ModifyPath"
+ Push "$INSTDIR\AddRemove.exe"
+ Call ConditionalAddToRegisty
+ !else
+ Push "NoModify"
+ Push "1"
+ Call ConditionalAddToRegisty
+ !endif
+
+ ; Optional registration
+ Push "DisplayIcon"
+ Push "$INSTDIR\@CPACK_NSIS_INSTALLED_ICON_NAME@"
+ Call ConditionalAddToRegisty
+ Push "HelpLink"
+ Push "@CPACK_NSIS_HELP_LINK@"
+ Call ConditionalAddToRegisty
+ Push "URLInfoAbout"
+ Push "@CPACK_NSIS_URL_INFO_ABOUT@"
+ Call ConditionalAddToRegisty
+ Push "Contact"
+ Push "@CPACK_NSIS_CONTACT@"
+ Call ConditionalAddToRegisty
+ !insertmacro MUI_INSTALLOPTIONS_READ $INSTALL_DESKTOP "NSIS.InstallOptions.ini" "Field 5" "State"
+ !insertmacro MUI_STARTMENU_WRITE_BEGIN Application
+
+ ;Create shortcuts
+ CreateDirectory "$SMPROGRAMS\$STARTMENU_FOLDER"
+ at CPACK_NSIS_CREATE_ICONS@
+ at CPACK_NSIS_CREATE_ICONS_EXTRA@
+ CreateShortCut "$SMPROGRAMS\$STARTMENU_FOLDER\Uninstall.lnk" "$INSTDIR\Uninstall.exe"
+
+ ; Write special uninstall registry entries
+ Push "StartMenu"
+ Push "$STARTMENU_FOLDER"
+ Call ConditionalAddToRegisty
+ Push "DoNotAddToPath"
+ Push "$DO_NOT_ADD_TO_PATH"
+ Call ConditionalAddToRegisty
+ Push "AddToPathAllUsers"
+ Push "$ADD_TO_PATH_ALL_USERS"
+ Call ConditionalAddToRegisty
+ Push "AddToPathCurrentUser"
+ Push "$ADD_TO_PATH_CURRENT_USER"
+ Call ConditionalAddToRegisty
+ Push "InstallToDesktop"
+ Push "$INSTALL_DESKTOP"
+ Call ConditionalAddToRegisty
+
+ at CPACK_NSIS_EXTRA_INSTALL_COMMANDS@
+
+ !insertmacro MUI_STARTMENU_WRITE_END
+
+SectionEnd
+
+Section "-Add to path"
+ Push $INSTDIR\bin
+ ;Read a value from an InstallOptions INI file
+ !insertmacro MUI_INSTALLOPTIONS_READ $DO_NOT_ADD_TO_PATH "NSIS.InstallOptions.ini" "Field 2" "State"
+ !insertmacro MUI_INSTALLOPTIONS_READ $ADD_TO_PATH_ALL_USERS "NSIS.InstallOptions.ini" "Field 3" "State"
+ !insertmacro MUI_INSTALLOPTIONS_READ $ADD_TO_PATH_CURRENT_USER "NSIS.InstallOptions.ini" "Field 4" "State"
+ StrCmp "@CPACK_NSIS_MODIFY_PATH@" "ON" 0 doNotAddToPath
+ StrCmp $DO_NOT_ADD_TO_PATH "1" doNotAddToPath 0
+ Call AddToPath
+ doNotAddToPath:
+SectionEnd
+
+;--------------------------------
+; Create custom pages
+Function InstallOptionsPage
+ !insertmacro MUI_HEADER_TEXT "Install Options" "Chose options for installing @CPACK_PACKAGE_INSTALL_DIRECTORY@"
+ !insertmacro MUI_INSTALLOPTIONS_DISPLAY "NSIS.InstallOptions.ini"
+
+FunctionEnd
+
+;--------------------------------
+; determine admin versus local install
+Function un.onInit
+
+ ClearErrors
+ UserInfo::GetName
+ IfErrors noLM
+ Pop $0
+ UserInfo::GetAccountType
+ Pop $1
+ StrCmp $1 "Admin" 0 +3
+ SetShellVarContext all
+ ;MessageBox MB_OK 'User "$0" is in the Admin group'
+ Goto done
+ StrCmp $1 "Power" 0 +3
+ SetShellVarContext all
+ ;MessageBox MB_OK 'User "$0" is in the Power Users group'
+ Goto done
+
+ noLM:
+ ;Get installation folder from registry if available
+
+ done:
+
+FunctionEnd
+
+;--- Add/Remove callback functions: ---
+!macro SectionList MacroName
+ ;This macro used to perform operation on multiple sections.
+ ;List all of your components in following manner here.
+ at CPACK_NSIS_COMPONENT_SECTION_LIST@
+!macroend
+
+Section -FinishComponents
+ ;Removes unselected components and writes component status to registry
+ !insertmacro SectionList "FinishSection"
+
+!ifdef CPACK_NSIS_ADD_REMOVE
+ ; Get the name of the installer executable
+ System::Call 'kernel32::GetModuleFileNameA(i 0, t .R0, i 1024) i r1'
+ StrCpy $R3 $R0
+
+ ; Strip off the last 13 characters, to see if we have AddRemove.exe
+ StrLen $R1 $R0
+ IntOp $R1 $R0 - 13
+ StrCpy $R2 $R0 13 $R1
+ StrCmp $R2 "AddRemove.exe" addremove_installed
+
+ ; We're not running AddRemove.exe, so install it
+ CopyFiles $R3 $INSTDIR\AddRemove.exe
+
+ addremove_installed:
+!endif
+SectionEnd
+;--- End of Add/Remove callback functions ---
+
+;--------------------------------
+; Component dependencies
+Function .onSelChange
+ !insertmacro SectionList MaybeSelectionChanged
+FunctionEnd
+
+;--------------------------------
+;Uninstaller Section
+
+Section "Uninstall"
+ ReadRegStr $START_MENU SHCTX \
+ "Software\Microsoft\Windows\CurrentVersion\Uninstall\@CPACK_PACKAGE_INSTALL_DIRECTORY@" "StartMenu"
+ ;MessageBox MB_OK "Start menu is in: $START_MENU"
+ ReadRegStr $DO_NOT_ADD_TO_PATH SHCTX \
+ "Software\Microsoft\Windows\CurrentVersion\Uninstall\@CPACK_PACKAGE_INSTALL_DIRECTORY@" "DoNotAddToPath"
+ ReadRegStr $ADD_TO_PATH_ALL_USERS SHCTX \
+ "Software\Microsoft\Windows\CurrentVersion\Uninstall\@CPACK_PACKAGE_INSTALL_DIRECTORY@" "AddToPathAllUsers"
+ ReadRegStr $ADD_TO_PATH_CURRENT_USER SHCTX \
+ "Software\Microsoft\Windows\CurrentVersion\Uninstall\@CPACK_PACKAGE_INSTALL_DIRECTORY@" "AddToPathCurrentUser"
+ ;MessageBox MB_OK "Add to path: $DO_NOT_ADD_TO_PATH all users: $ADD_TO_PATH_ALL_USERS"
+ ReadRegStr $INSTALL_DESKTOP SHCTX \
+ "Software\Microsoft\Windows\CurrentVersion\Uninstall\@CPACK_PACKAGE_INSTALL_DIRECTORY@" "InstallToDesktop"
+ ;MessageBox MB_OK "Install to desktop: $INSTALL_DESKTOP "
+
+ at CPACK_NSIS_EXTRA_UNINSTALL_COMMANDS@
+
+ ;Remove files we installed.
+ ;Keep the list of directories here in sync with the File commands above.
+ at CPACK_NSIS_DELETE_FILES@
+ at CPACK_NSIS_DELETE_DIRECTORIES@
+
+!ifdef CPACK_NSIS_ADD_REMOVE
+ ;Remove the add/remove program
+ Delete "$INSTDIR\AddRemove.exe"
+!endif
+
+ ;Remove the uninstaller itself.
+ Delete "$INSTDIR\Uninstall.exe"
+ DeleteRegKey SHCTX "Software\Microsoft\Windows\CurrentVersion\Uninstall\@CPACK_PACKAGE_INSTALL_DIRECTORY@"
+
+ ;Remove the installation directory if it is empty.
+ RMDir "$INSTDIR"
+
+ ; Remove the registry entries.
+ DeleteRegKey SHCTX "Software\@CPACK_PACKAGE_VENDOR@\@CPACK_PACKAGE_INSTALL_REGISTRY_KEY@"
+
+ ; Removes all optional components
+ !insertmacro SectionList "RemoveSection"
+
+ !insertmacro MUI_STARTMENU_GETFOLDER Application $MUI_TEMP
+
+ Delete "$SMPROGRAMS\$MUI_TEMP\Uninstall.lnk"
+ at CPACK_NSIS_DELETE_ICONS@
+ at CPACK_NSIS_DELETE_ICONS_EXTRA@
+
+ ;Delete empty start menu parent diretories
+ StrCpy $MUI_TEMP "$SMPROGRAMS\$MUI_TEMP"
+
+ startMenuDeleteLoop:
+ ClearErrors
+ RMDir $MUI_TEMP
+ GetFullPathName $MUI_TEMP "$MUI_TEMP\.."
+
+ IfErrors startMenuDeleteLoopDone
+
+ StrCmp "$MUI_TEMP" "$SMPROGRAMS" startMenuDeleteLoopDone startMenuDeleteLoop
+ startMenuDeleteLoopDone:
+
+ ; If the user changed the shortcut, then untinstall may not work. This should
+ ; try to fix it.
+ StrCpy $MUI_TEMP "$START_MENU"
+ Delete "$SMPROGRAMS\$MUI_TEMP\Uninstall.lnk"
+ at CPACK_NSIS_DELETE_ICONS_EXTRA@
+
+ ;Delete empty start menu parent diretories
+ StrCpy $MUI_TEMP "$SMPROGRAMS\$MUI_TEMP"
+
+ secondStartMenuDeleteLoop:
+ ClearErrors
+ RMDir $MUI_TEMP
+ GetFullPathName $MUI_TEMP "$MUI_TEMP\.."
+
+ IfErrors secondStartMenuDeleteLoopDone
+
+ StrCmp "$MUI_TEMP" "$SMPROGRAMS" secondStartMenuDeleteLoopDone secondStartMenuDeleteLoop
+ secondStartMenuDeleteLoopDone:
+
+ DeleteRegKey /ifempty SHCTX "Software\@CPACK_PACKAGE_VENDOR@\@CPACK_PACKAGE_INSTALL_REGISTRY_KEY@"
+
+ Push $INSTDIR\bin
+ StrCmp $DO_NOT_ADD_TO_PATH_ "1" doNotRemoveFromPath 0
+ Call un.RemoveFromPath
+ doNotRemoveFromPath:
+SectionEnd
+
+;--------------------------------
+; determine admin versus local install
+; Is install for "AllUsers" or "JustMe"?
+; Default to "JustMe" - set to "AllUsers" if admin or on Win9x
+; This function is used for the very first "custom page" of the installer.
+; This custom page does not show up visibly, but it executes prior to the
+; first visible page and sets up $INSTDIR properly...
+; Choose different default installation folder based on SV_ALLUSERS...
+; "Program Files" for AllUsers, "My Documents" for JustMe...
+
+Function .onInit
+ ; Reads components status for registry
+ !insertmacro SectionList "InitSection"
+
+ StrCpy $SV_ALLUSERS "JustMe"
+ StrCpy $INSTDIR "$DOCUMENTS\@CPACK_PACKAGE_INSTALL_DIRECTORY@"
+
+ ClearErrors
+ UserInfo::GetName
+ IfErrors noLM
+ Pop $0
+ UserInfo::GetAccountType
+ Pop $1
+ StrCmp $1 "Admin" 0 +3
+ SetShellVarContext all
+ ;MessageBox MB_OK 'User "$0" is in the Admin group'
+ StrCpy $SV_ALLUSERS "AllUsers"
+ Goto done
+ StrCmp $1 "Power" 0 +3
+ SetShellVarContext all
+ ;MessageBox MB_OK 'User "$0" is in the Power Users group'
+ StrCpy $SV_ALLUSERS "AllUsers"
+ Goto done
+
+ noLM:
+ StrCpy $SV_ALLUSERS "AllUsers"
+ ;Get installation folder from registry if available
+
+ done:
+ StrCmp $SV_ALLUSERS "AllUsers" 0 +2
+ StrCpy $INSTDIR "$PROGRAMFILES\@CPACK_PACKAGE_INSTALL_DIRECTORY@"
+
+ StrCmp "@CPACK_NSIS_MODIFY_PATH@" "ON" 0 noOptionsPage
+ !insertmacro MUI_INSTALLOPTIONS_EXTRACT "NSIS.InstallOptions.ini"
+
+ noOptionsPage:
+FunctionEnd
diff --git a/copying.txt b/copying.txt
new file mode 100644
index 0000000..10926e8
--- /dev/null
+++ b/copying.txt
@@ -0,0 +1,675 @@
+ GNU GENERAL PUBLIC LICENSE
+ Version 3, 29 June 2007
+
+ Copyright (C) 2007 Free Software Foundation, Inc. <http://fsf.org/>
+ Everyone is permitted to copy and distribute verbatim copies
+ of this license document, but changing it is not allowed.
+
+ Preamble
+
+ The GNU General Public License is a free, copyleft license for
+software and other kinds of works.
+
+ The licenses for most software and other practical works are designed
+to take away your freedom to share and change the works. By contrast,
+the GNU General Public License is intended to guarantee your freedom to
+share and change all versions of a program--to make sure it remains free
+software for all its users. We, the Free Software Foundation, use the
+GNU General Public License for most of our software; it applies also to
+any other work released this way by its authors. You can apply it to
+your programs, too.
+
+ When we speak of free software, we are referring to freedom, not
+price. Our General Public Licenses are designed to make sure that you
+have the freedom to distribute copies of free software (and charge for
+them if you wish), that you receive source code or can get it if you
+want it, that you can change the software or use pieces of it in new
+free programs, and that you know you can do these things.
+
+ To protect your rights, we need to prevent others from denying you
+these rights or asking you to surrender the rights. Therefore, you have
+certain responsibilities if you distribute copies of the software, or if
+you modify it: responsibilities to respect the freedom of others.
+
+ For example, if you distribute copies of such a program, whether
+gratis or for a fee, you must pass on to the recipients the same
+freedoms that you received. You must make sure that they, too, receive
+or can get the source code. And you must show them these terms so they
+know their rights.
+
+ Developers that use the GNU GPL protect your rights with two steps:
+(1) assert copyright on the software, and (2) offer you this License
+giving you legal permission to copy, distribute and/or modify it.
+
+ For the developers' and authors' protection, the GPL clearly explains
+that there is no warranty for this free software. For both users' and
+authors' sake, the GPL requires that modified versions be marked as
+changed, so that their problems will not be attributed erroneously to
+authors of previous versions.
+
+ Some devices are designed to deny users access to install or run
+modified versions of the software inside them, although the manufacturer
+can do so. This is fundamentally incompatible with the aim of
+protecting users' freedom to change the software. The systematic
+pattern of such abuse occurs in the area of products for individuals to
+use, which is precisely where it is most unacceptable. Therefore, we
+have designed this version of the GPL to prohibit the practice for those
+products. If such problems arise substantially in other domains, we
+stand ready to extend this provision to those domains in future versions
+of the GPL, as needed to protect the freedom of users.
+
+ Finally, every program is threatened constantly by software patents.
+States should not allow patents to restrict development and use of
+software on general-purpose computers, but in those that do, we wish to
+avoid the special danger that patents applied to a free program could
+make it effectively proprietary. To prevent this, the GPL assures that
+patents cannot be used to render the program non-free.
+
+ The precise terms and conditions for copying, distribution and
+modification follow.
+
+ TERMS AND CONDITIONS
+
+ 0. Definitions.
+
+ "This License" refers to version 3 of the GNU General Public License.
+
+ "Copyright" also means copyright-like laws that apply to other kinds of
+works, such as semiconductor masks.
+
+ "The Program" refers to any copyrightable work licensed under this
+License. Each licensee is addressed as "you". "Licensees" and
+"recipients" may be individuals or organizations.
+
+ To "modify" a work means to copy from or adapt all or part of the work
+in a fashion requiring copyright permission, other than the making of an
+exact copy. The resulting work is called a "modified version" of the
+earlier work or a work "based on" the earlier work.
+
+ A "covered work" means either the unmodified Program or a work based
+on the Program.
+
+ To "propagate" a work means to do anything with it that, without
+permission, would make you directly or secondarily liable for
+infringement under applicable copyright law, except executing it on a
+computer or modifying a private copy. Propagation includes copying,
+distribution (with or without modification), making available to the
+public, and in some countries other activities as well.
+
+ To "convey" a work means any kind of propagation that enables other
+parties to make or receive copies. Mere interaction with a user through
+a computer network, with no transfer of a copy, is not conveying.
+
+ An interactive user interface displays "Appropriate Legal Notices"
+to the extent that it includes a convenient and prominently visible
+feature that (1) displays an appropriate copyright notice, and (2)
+tells the user that there is no warranty for the work (except to the
+extent that warranties are provided), that licensees may convey the
+work under this License, and how to view a copy of this License. If
+the interface presents a list of user commands or options, such as a
+menu, a prominent item in the list meets this criterion.
+
+ 1. Source Code.
+
+ The "source code" for a work means the preferred form of the work
+for making modifications to it. "Object code" means any non-source
+form of a work.
+
+ A "Standard Interface" means an interface that either is an official
+standard defined by a recognized standards body, or, in the case of
+interfaces specified for a particular programming language, one that
+is widely used among developers working in that language.
+
+ The "System Libraries" of an executable work include anything, other
+than the work as a whole, that (a) is included in the normal form of
+packaging a Major Component, but which is not part of that Major
+Component, and (b) serves only to enable use of the work with that
+Major Component, or to implement a Standard Interface for which an
+implementation is available to the public in source code form. A
+"Major Component", in this context, means a major essential component
+(kernel, window system, and so on) of the specific operating system
+(if any) on which the executable work runs, or a compiler used to
+produce the work, or an object code interpreter used to run it.
+
+ The "Corresponding Source" for a work in object code form means all
+the source code needed to generate, install, and (for an executable
+work) run the object code and to modify the work, including scripts to
+control those activities. However, it does not include the work's
+System Libraries, or general-purpose tools or generally available free
+programs which are used unmodified in performing those activities but
+which are not part of the work. For example, Corresponding Source
+includes interface definition files associated with source files for
+the work, and the source code for shared libraries and dynamically
+linked subprograms that the work is specifically designed to require,
+such as by intimate data communication or control flow between those
+subprograms and other parts of the work.
+
+ The Corresponding Source need not include anything that users
+can regenerate automatically from other parts of the Corresponding
+Source.
+
+ The Corresponding Source for a work in source code form is that
+same work.
+
+ 2. Basic Permissions.
+
+ All rights granted under this License are granted for the term of
+copyright on the Program, and are irrevocable provided the stated
+conditions are met. This License explicitly affirms your unlimited
+permission to run the unmodified Program. The output from running a
+covered work is covered by this License only if the output, given its
+content, constitutes a covered work. This License acknowledges your
+rights of fair use or other equivalent, as provided by copyright law.
+
+ You may make, run and propagate covered works that you do not
+convey, without conditions so long as your license otherwise remains
+in force. You may convey covered works to others for the sole purpose
+of having them make modifications exclusively for you, or provide you
+with facilities for running those works, provided that you comply with
+the terms of this License in conveying all material for which you do
+not control copyright. Those thus making or running the covered works
+for you must do so exclusively on your behalf, under your direction
+and control, on terms that prohibit them from making any copies of
+your copyrighted material outside their relationship with you.
+
+ Conveying under any other circumstances is permitted solely under
+the conditions stated below. Sublicensing is not allowed; section 10
+makes it unnecessary.
+
+ 3. Protecting Users' Legal Rights From Anti-Circumvention Law.
+
+ No covered work shall be deemed part of an effective technological
+measure under any applicable law fulfilling obligations under article
+11 of the WIPO copyright treaty adopted on 20 December 1996, or
+similar laws prohibiting or restricting circumvention of such
+measures.
+
+ When you convey a covered work, you waive any legal power to forbid
+circumvention of technological measures to the extent such circumvention
+is effected by exercising rights under this License with respect to
+the covered work, and you disclaim any intention to limit operation or
+modification of the work as a means of enforcing, against the work's
+users, your or third parties' legal rights to forbid circumvention of
+technological measures.
+
+ 4. Conveying Verbatim Copies.
+
+ You may convey verbatim copies of the Program's source code as you
+receive it, in any medium, provided that you conspicuously and
+appropriately publish on each copy an appropriate copyright notice;
+keep intact all notices stating that this License and any
+non-permissive terms added in accord with section 7 apply to the code;
+keep intact all notices of the absence of any warranty; and give all
+recipients a copy of this License along with the Program.
+
+ You may charge any price or no price for each copy that you convey,
+and you may offer support or warranty protection for a fee.
+
+ 5. Conveying Modified Source Versions.
+
+ You may convey a work based on the Program, or the modifications to
+produce it from the Program, in the form of source code under the
+terms of section 4, provided that you also meet all of these conditions:
+
+ a) The work must carry prominent notices stating that you modified
+ it, and giving a relevant date.
+
+ b) The work must carry prominent notices stating that it is
+ released under this License and any conditions added under section
+ 7. This requirement modifies the requirement in section 4 to
+ "keep intact all notices".
+
+ c) You must license the entire work, as a whole, under this
+ License to anyone who comes into possession of a copy. This
+ License will therefore apply, along with any applicable section 7
+ additional terms, to the whole of the work, and all its parts,
+ regardless of how they are packaged. This License gives no
+ permission to license the work in any other way, but it does not
+ invalidate such permission if you have separately received it.
+
+ d) If the work has interactive user interfaces, each must display
+ Appropriate Legal Notices; however, if the Program has interactive
+ interfaces that do not display Appropriate Legal Notices, your
+ work need not make them do so.
+
+ A compilation of a covered work with other separate and independent
+works, which are not by their nature extensions of the covered work,
+and which are not combined with it such as to form a larger program,
+in or on a volume of a storage or distribution medium, is called an
+"aggregate" if the compilation and its resulting copyright are not
+used to limit the access or legal rights of the compilation's users
+beyond what the individual works permit. Inclusion of a covered work
+in an aggregate does not cause this License to apply to the other
+parts of the aggregate.
+
+ 6. Conveying Non-Source Forms.
+
+ You may convey a covered work in object code form under the terms
+of sections 4 and 5, provided that you also convey the
+machine-readable Corresponding Source under the terms of this License,
+in one of these ways:
+
+ a) Convey the object code in, or embodied in, a physical product
+ (including a physical distribution medium), accompanied by the
+ Corresponding Source fixed on a durable physical medium
+ customarily used for software interchange.
+
+ b) Convey the object code in, or embodied in, a physical product
+ (including a physical distribution medium), accompanied by a
+ written offer, valid for at least three years and valid for as
+ long as you offer spare parts or customer support for that product
+ model, to give anyone who possesses the object code either (1) a
+ copy of the Corresponding Source for all the software in the
+ product that is covered by this License, on a durable physical
+ medium customarily used for software interchange, for a price no
+ more than your reasonable cost of physically performing this
+ conveying of source, or (2) access to copy the
+ Corresponding Source from a network server at no charge.
+
+ c) Convey individual copies of the object code with a copy of the
+ written offer to provide the Corresponding Source. This
+ alternative is allowed only occasionally and noncommercially, and
+ only if you received the object code with such an offer, in accord
+ with subsection 6b.
+
+ d) Convey the object code by offering access from a designated
+ place (gratis or for a charge), and offer equivalent access to the
+ Corresponding Source in the same way through the same place at no
+ further charge. You need not require recipients to copy the
+ Corresponding Source along with the object code. If the place to
+ copy the object code is a network server, the Corresponding Source
+ may be on a different server (operated by you or a third party)
+ that supports equivalent copying facilities, provided you maintain
+ clear directions next to the object code saying where to find the
+ Corresponding Source. Regardless of what server hosts the
+ Corresponding Source, you remain obligated to ensure that it is
+ available for as long as needed to satisfy these requirements.
+
+ e) Convey the object code using peer-to-peer transmission, provided
+ you inform other peers where the object code and Corresponding
+ Source of the work are being offered to the general public at no
+ charge under subsection 6d.
+
+ A separable portion of the object code, whose source code is excluded
+from the Corresponding Source as a System Library, need not be
+included in conveying the object code work.
+
+ A "User Product" is either (1) a "consumer product", which means any
+tangible personal property which is normally used for personal, family,
+or household purposes, or (2) anything designed or sold for incorporation
+into a dwelling. In determining whether a product is a consumer product,
+doubtful cases shall be resolved in favor of coverage. For a particular
+product received by a particular user, "normally used" refers to a
+typical or common use of that class of product, regardless of the status
+of the particular user or of the way in which the particular user
+actually uses, or expects or is expected to use, the product. A product
+is a consumer product regardless of whether the product has substantial
+commercial, industrial or non-consumer uses, unless such uses represent
+the only significant mode of use of the product.
+
+ "Installation Information" for a User Product means any methods,
+procedures, authorization keys, or other information required to install
+and execute modified versions of a covered work in that User Product from
+a modified version of its Corresponding Source. The information must
+suffice to ensure that the continued functioning of the modified object
+code is in no case prevented or interfered with solely because
+modification has been made.
+
+ If you convey an object code work under this section in, or with, or
+specifically for use in, a User Product, and the conveying occurs as
+part of a transaction in which the right of possession and use of the
+User Product is transferred to the recipient in perpetuity or for a
+fixed term (regardless of how the transaction is characterized), the
+Corresponding Source conveyed under this section must be accompanied
+by the Installation Information. But this requirement does not apply
+if neither you nor any third party retains the ability to install
+modified object code on the User Product (for example, the work has
+been installed in ROM).
+
+ The requirement to provide Installation Information does not include a
+requirement to continue to provide support service, warranty, or updates
+for a work that has been modified or installed by the recipient, or for
+the User Product in which it has been modified or installed. Access to a
+network may be denied when the modification itself materially and
+adversely affects the operation of the network or violates the rules and
+protocols for communication across the network.
+
+ Corresponding Source conveyed, and Installation Information provided,
+in accord with this section must be in a format that is publicly
+documented (and with an implementation available to the public in
+source code form), and must require no special password or key for
+unpacking, reading or copying.
+
+ 7. Additional Terms.
+
+ "Additional permissions" are terms that supplement the terms of this
+License by making exceptions from one or more of its conditions.
+Additional permissions that are applicable to the entire Program shall
+be treated as though they were included in this License, to the extent
+that they are valid under applicable law. If additional permissions
+apply only to part of the Program, that part may be used separately
+under those permissions, but the entire Program remains governed by
+this License without regard to the additional permissions.
+
+ When you convey a copy of a covered work, you may at your option
+remove any additional permissions from that copy, or from any part of
+it. (Additional permissions may be written to require their own
+removal in certain cases when you modify the work.) You may place
+additional permissions on material, added by you to a covered work,
+for which you have or can give appropriate copyright permission.
+
+ Notwithstanding any other provision of this License, for material you
+add to a covered work, you may (if authorized by the copyright holders of
+that material) supplement the terms of this License with terms:
+
+ a) Disclaiming warranty or limiting liability differently from the
+ terms of sections 15 and 16 of this License; or
+
+ b) Requiring preservation of specified reasonable legal notices or
+ author attributions in that material or in the Appropriate Legal
+ Notices displayed by works containing it; or
+
+ c) Prohibiting misrepresentation of the origin of that material, or
+ requiring that modified versions of such material be marked in
+ reasonable ways as different from the original version; or
+
+ d) Limiting the use for publicity purposes of names of licensors or
+ authors of the material; or
+
+ e) Declining to grant rights under trademark law for use of some
+ trade names, trademarks, or service marks; or
+
+ f) Requiring indemnification of licensors and authors of that
+ material by anyone who conveys the material (or modified versions of
+ it) with contractual assumptions of liability to the recipient, for
+ any liability that these contractual assumptions directly impose on
+ those licensors and authors.
+
+ All other non-permissive additional terms are considered "further
+restrictions" within the meaning of section 10. If the Program as you
+received it, or any part of it, contains a notice stating that it is
+governed by this License along with a term that is a further
+restriction, you may remove that term. If a license document contains
+a further restriction but permits relicensing or conveying under this
+License, you may add to a covered work material governed by the terms
+of that license document, provided that the further restriction does
+not survive such relicensing or conveying.
+
+ If you add terms to a covered work in accord with this section, you
+must place, in the relevant source files, a statement of the
+additional terms that apply to those files, or a notice indicating
+where to find the applicable terms.
+
+ Additional terms, permissive or non-permissive, may be stated in the
+form of a separately written license, or stated as exceptions;
+the above requirements apply either way.
+
+ 8. Termination.
+
+ You may not propagate or modify a covered work except as expressly
+provided under this License. Any attempt otherwise to propagate or
+modify it is void, and will automatically terminate your rights under
+this License (including any patent licenses granted under the third
+paragraph of section 11).
+
+ However, if you cease all violation of this License, then your
+license from a particular copyright holder is reinstated (a)
+provisionally, unless and until the copyright holder explicitly and
+finally terminates your license, and (b) permanently, if the copyright
+holder fails to notify you of the violation by some reasonable means
+prior to 60 days after the cessation.
+
+ Moreover, your license from a particular copyright holder is
+reinstated permanently if the copyright holder notifies you of the
+violation by some reasonable means, this is the first time you have
+received notice of violation of this License (for any work) from that
+copyright holder, and you cure the violation prior to 30 days after
+your receipt of the notice.
+
+ Termination of your rights under this section does not terminate the
+licenses of parties who have received copies or rights from you under
+this License. If your rights have been terminated and not permanently
+reinstated, you do not qualify to receive new licenses for the same
+material under section 10.
+
+ 9. Acceptance Not Required for Having Copies.
+
+ You are not required to accept this License in order to receive or
+run a copy of the Program. Ancillary propagation of a covered work
+occurring solely as a consequence of using peer-to-peer transmission
+to receive a copy likewise does not require acceptance. However,
+nothing other than this License grants you permission to propagate or
+modify any covered work. These actions infringe copyright if you do
+not accept this License. Therefore, by modifying or propagating a
+covered work, you indicate your acceptance of this License to do so.
+
+ 10. Automatic Licensing of Downstream Recipients.
+
+ Each time you convey a covered work, the recipient automatically
+receives a license from the original licensors, to run, modify and
+propagate that work, subject to this License. You are not responsible
+for enforcing compliance by third parties with this License.
+
+ An "entity transaction" is a transaction transferring control of an
+organization, or substantially all assets of one, or subdividing an
+organization, or merging organizations. If propagation of a covered
+work results from an entity transaction, each party to that
+transaction who receives a copy of the work also receives whatever
+licenses to the work the party's predecessor in interest had or could
+give under the previous paragraph, plus a right to possession of the
+Corresponding Source of the work from the predecessor in interest, if
+the predecessor has it or can get it with reasonable efforts.
+
+ You may not impose any further restrictions on the exercise of the
+rights granted or affirmed under this License. For example, you may
+not impose a license fee, royalty, or other charge for exercise of
+rights granted under this License, and you may not initiate litigation
+(including a cross-claim or counterclaim in a lawsuit) alleging that
+any patent claim is infringed by making, using, selling, offering for
+sale, or importing the Program or any portion of it.
+
+ 11. Patents.
+
+ A "contributor" is a copyright holder who authorizes use under this
+License of the Program or a work on which the Program is based. The
+work thus licensed is called the contributor's "contributor version".
+
+ A contributor's "essential patent claims" are all patent claims
+owned or controlled by the contributor, whether already acquired or
+hereafter acquired, that would be infringed by some manner, permitted
+by this License, of making, using, or selling its contributor version,
+but do not include claims that would be infringed only as a
+consequence of further modification of the contributor version. For
+purposes of this definition, "control" includes the right to grant
+patent sublicenses in a manner consistent with the requirements of
+this License.
+
+ Each contributor grants you a non-exclusive, worldwide, royalty-free
+patent license under the contributor's essential patent claims, to
+make, use, sell, offer for sale, import and otherwise run, modify and
+propagate the contents of its contributor version.
+
+ In the following three paragraphs, a "patent license" is any express
+agreement or commitment, however denominated, not to enforce a patent
+(such as an express permission to practice a patent or covenant not to
+sue for patent infringement). To "grant" such a patent license to a
+party means to make such an agreement or commitment not to enforce a
+patent against the party.
+
+ If you convey a covered work, knowingly relying on a patent license,
+and the Corresponding Source of the work is not available for anyone
+to copy, free of charge and under the terms of this License, through a
+publicly available network server or other readily accessible means,
+then you must either (1) cause the Corresponding Source to be so
+available, or (2) arrange to deprive yourself of the benefit of the
+patent license for this particular work, or (3) arrange, in a manner
+consistent with the requirements of this License, to extend the patent
+license to downstream recipients. "Knowingly relying" means you have
+actual knowledge that, but for the patent license, your conveying the
+covered work in a country, or your recipient's use of the covered work
+in a country, would infringe one or more identifiable patents in that
+country that you have reason to believe are valid.
+
+ If, pursuant to or in connection with a single transaction or
+arrangement, you convey, or propagate by procuring conveyance of, a
+covered work, and grant a patent license to some of the parties
+receiving the covered work authorizing them to use, propagate, modify
+or convey a specific copy of the covered work, then the patent license
+you grant is automatically extended to all recipients of the covered
+work and works based on it.
+
+ A patent license is "discriminatory" if it does not include within
+the scope of its coverage, prohibits the exercise of, or is
+conditioned on the non-exercise of one or more of the rights that are
+specifically granted under this License. You may not convey a covered
+work if you are a party to an arrangement with a third party that is
+in the business of distributing software, under which you make payment
+to the third party based on the extent of your activity of conveying
+the work, and under which the third party grants, to any of the
+parties who would receive the covered work from you, a discriminatory
+patent license (a) in connection with copies of the covered work
+conveyed by you (or copies made from those copies), or (b) primarily
+for and in connection with specific products or compilations that
+contain the covered work, unless you entered into that arrangement,
+or that patent license was granted, prior to 28 March 2007.
+
+ Nothing in this License shall be construed as excluding or limiting
+any implied license or other defenses to infringement that may
+otherwise be available to you under applicable patent law.
+
+ 12. No Surrender of Others' Freedom.
+
+ If conditions are imposed on you (whether by court order, agreement or
+otherwise) that contradict the conditions of this License, they do not
+excuse you from the conditions of this License. If you cannot convey a
+covered work so as to satisfy simultaneously your obligations under this
+License and any other pertinent obligations, then as a consequence you may
+not convey it at all. For example, if you agree to terms that obligate you
+to collect a royalty for further conveying from those to whom you convey
+the Program, the only way you could satisfy both those terms and this
+License would be to refrain entirely from conveying the Program.
+
+ 13. Use with the GNU Affero General Public License.
+
+ Notwithstanding any other provision of this License, you have
+permission to link or combine any covered work with a work licensed
+under version 3 of the GNU Affero General Public License into a single
+combined work, and to convey the resulting work. The terms of this
+License will continue to apply to the part which is the covered work,
+but the special requirements of the GNU Affero General Public License,
+section 13, concerning interaction through a network will apply to the
+combination as such.
+
+ 14. Revised Versions of this License.
+
+ The Free Software Foundation may publish revised and/or new versions of
+the GNU General Public License from time to time. Such new versions will
+be similar in spirit to the present version, but may differ in detail to
+address new problems or concerns.
+
+ Each version is given a distinguishing version number. If the
+Program specifies that a certain numbered version of the GNU General
+Public License "or any later version" applies to it, you have the
+option of following the terms and conditions either of that numbered
+version or of any later version published by the Free Software
+Foundation. If the Program does not specify a version number of the
+GNU General Public License, you may choose any version ever published
+by the Free Software Foundation.
+
+ If the Program specifies that a proxy can decide which future
+versions of the GNU General Public License can be used, that proxy's
+public statement of acceptance of a version permanently authorizes you
+to choose that version for the Program.
+
+ Later license versions may give you additional or different
+permissions. However, no additional obligations are imposed on any
+author or copyright holder as a result of your choosing to follow a
+later version.
+
+ 15. Disclaimer of Warranty.
+
+ THERE IS NO WARRANTY FOR THE PROGRAM, TO THE EXTENT PERMITTED BY
+APPLICABLE LAW. EXCEPT WHEN OTHERWISE STATED IN WRITING THE COPYRIGHT
+HOLDERS AND/OR OTHER PARTIES PROVIDE THE PROGRAM "AS IS" WITHOUT WARRANTY
+OF ANY KIND, EITHER EXPRESSED OR IMPLIED, INCLUDING, BUT NOT LIMITED TO,
+THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
+PURPOSE. THE ENTIRE RISK AS TO THE QUALITY AND PERFORMANCE OF THE PROGRAM
+IS WITH YOU. SHOULD THE PROGRAM PROVE DEFECTIVE, YOU ASSUME THE COST OF
+ALL NECESSARY SERVICING, REPAIR OR CORRECTION.
+
+ 16. Limitation of Liability.
+
+ IN NO EVENT UNLESS REQUIRED BY APPLICABLE LAW OR AGREED TO IN WRITING
+WILL ANY COPYRIGHT HOLDER, OR ANY OTHER PARTY WHO MODIFIES AND/OR CONVEYS
+THE PROGRAM AS PERMITTED ABOVE, BE LIABLE TO YOU FOR DAMAGES, INCLUDING ANY
+GENERAL, SPECIAL, INCIDENTAL OR CONSEQUENTIAL DAMAGES ARISING OUT OF THE
+USE OR INABILITY TO USE THE PROGRAM (INCLUDING BUT NOT LIMITED TO LOSS OF
+DATA OR DATA BEING RENDERED INACCURATE OR LOSSES SUSTAINED BY YOU OR THIRD
+PARTIES OR A FAILURE OF THE PROGRAM TO OPERATE WITH ANY OTHER PROGRAMS),
+EVEN IF SUCH HOLDER OR OTHER PARTY HAS BEEN ADVISED OF THE POSSIBILITY OF
+SUCH DAMAGES.
+
+ 17. Interpretation of Sections 15 and 16.
+
+ If the disclaimer of warranty and limitation of liability provided
+above cannot be given local legal effect according to their terms,
+reviewing courts shall apply local law that most closely approximates
+an absolute waiver of all civil liability in connection with the
+Program, unless a warranty or assumption of liability accompanies a
+copy of the Program in return for a fee.
+
+ END OF TERMS AND CONDITIONS
+
+ How to Apply These Terms to Your New Programs
+
+ If you develop a new program, and you want it to be of the greatest
+possible use to the public, the best way to achieve this is to make it
+free software which everyone can redistribute and change under these terms.
+
+ To do so, attach the following notices to the program. It is safest
+to attach them to the start of each source file to most effectively
+state the exclusion of warranty; and each file should have at least
+the "copyright" line and a pointer to where the full notice is found.
+
+ <one line to give the program's name and a brief idea of what it does.>
+ Copyright (C) <year> <name of author>
+
+ This program is free software: you can redistribute it and/or modify
+ it under the terms of the GNU General Public License as published by
+ the Free Software Foundation, either version 3 of the License, or
+ (at your option) any later version.
+
+ This program is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ GNU General Public License for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with this program. If not, see <http://www.gnu.org/licenses/>.
+
+Also add information on how to contact you by electronic and paper mail.
+
+ If the program does terminal interaction, make it output a short
+notice like this when it starts in an interactive mode:
+
+ <program> Copyright (C) <year> <name of author>
+ This program comes with ABSOLUTELY NO WARRANTY; for details type `show w'.
+ This is free software, and you are welcome to redistribute it
+ under certain conditions; type `show c' for details.
+
+The hypothetical commands `show w' and `show c' should show the appropriate
+parts of the General Public License. Of course, your program's commands
+might be different; for a GUI interface, you would use an "about box".
+
+ You should also get your employer (if you work as a programmer) or school,
+if any, to sign a "copyright disclaimer" for the program, if necessary.
+For more information on this, and how to apply and follow the GNU GPL, see
+<http://www.gnu.org/licenses/>.
+
+ The GNU General Public License does not permit incorporating your program
+into proprietary programs. If your program is a subroutine library, you
+may consider it more useful to permit linking proprietary applications with
+the library. If this is what you want to do, use the GNU Lesser General
+Public License instead of this License. But first, please read
+<http://www.gnu.org/philosophy/why-not-lgpl.html>.
+
diff --git a/debian/changelog b/debian/changelog
deleted file mode 100644
index f7554f2..0000000
--- a/debian/changelog
+++ /dev/null
@@ -1,5 +0,0 @@
-ngila (1.3-1) UNRELEASED; urgency=low
-
- * Initial release (Closes: #439996)
-
- -- Andreas Tille <tille at debian.org> Tue, 04 Feb 2014 14:36:01 +0100
diff --git a/debian/compat b/debian/compat
deleted file mode 100644
index ec63514..0000000
--- a/debian/compat
+++ /dev/null
@@ -1 +0,0 @@
-9
diff --git a/debian/control b/debian/control
deleted file mode 100644
index ff36481..0000000
--- a/debian/control
+++ /dev/null
@@ -1,32 +0,0 @@
-Source: ngila
-Maintainer: Debian Med Packaging Team <debian-med-packaging at lists.alioth.debian.org>
-Uploaders: Andreas Tille <tille at debian.org>,
- Daniel Barker <db60 at st-andrews.ac.uk>
-Section: science
-Priority: optional
-Build-Depends: debhelper (>= 9)
-Standards-Version: 3.9.5
-Vcs-Browser: http://anonscm.debian.org/viewvc/debian-med/trunk/packages/ngila/trunk/
-Vcs-Svn: svn://anonscm.debian.org/debian-med/trunk/packages/ngila/trunk/
-Homepage: http://scit.us/projects/ngila/
-
-Package: ngila
-Architecture: any
-Depends: ${shlibs:Depends},
- ${misc:Depends}
-Description: global pairwise alignments with logarithmic and affine gap costs
- Ngila is an application that will find the best alignment of a pair
- of sequences using log-affine gap costs, which are the most
- biologically realistic gap costs.
- .
- Ngila implements the Miller and Myers (1988) algorithm in order to
- find a least costly global alignment of two sequences given homology
- costs and a gap cost. Two versions of the algorithm are
- included: holistic and divide-and-conquer. The former is faster but
- the latter utilizes less memory. Ngila starts with the
- divide-and-conquer method but switches to the holistic method for
- subsequences smaller than a user-established threshold. This improves
- its speed without substantially increasing memory requirements. Ngila
- also allows users to assign costs to end gaps that are smaller than
- costs for internal gaps. This is important for aligning using the
- free-end-gap method.
diff --git a/debian/copyright b/debian/copyright
deleted file mode 100644
index 4bd03f2..0000000
--- a/debian/copyright
+++ /dev/null
@@ -1,10 +0,0 @@
-Format: http://www.debian.org/doc/packaging-manuals/copyright-format/1.0/
-Source: http://scit.us/projects/ngila/downloads
-
-Files: *
-Copyright: © 2008-2011 Reed A. Cartwright <racartwright at rice.edu>
-License: GPLv3
-
-Files: debian/*
-Copyright: © 2012 Andreas Tille <tille at debian.org>
-License: GPLv3
diff --git a/debian/rules b/debian/rules
deleted file mode 100755
index ac62313..0000000
--- a/debian/rules
+++ /dev/null
@@ -1,10 +0,0 @@
-#!/usr/bin/make -f
-
-# DH_VERBOSE := 1
-
-%:
- dh $@
-
-get-orig-source:
- # Uostream does not provide reasonable versioning :-(
- wget http://scit.us/projects/files/ngila/Releases/ngila-release.tar.gz
diff --git a/debian/source/format b/debian/source/format
deleted file mode 100644
index 163aaf8..0000000
--- a/debian/source/format
+++ /dev/null
@@ -1 +0,0 @@
-3.0 (quilt)
diff --git a/debian/upstream/metadata b/debian/upstream/metadata
deleted file mode 100644
index ab97e4a..0000000
--- a/debian/upstream/metadata
+++ /dev/null
@@ -1,12 +0,0 @@
-Reference:
- Author: Reed A. Cartwright
- Title: "Ngila: global pairwise alignments with logarithmic and affine gap costs"
- Journal: Bioinformatics
- Year: 2007
- Volume: 23
- Number: 11
- Pages: 1427-1428
- DOI: 10.1093/bioinformatics/btm095
- PMID: 17387111
- URL: http://bioinformatics.oxfordjournals.org/content/23/11/1427
- eprint: http://bioinformatics.oxfordjournals.org/content/23/11/1427.full.pdf+html
diff --git a/debian/watch b/debian/watch
deleted file mode 100644
index 1390b91..0000000
--- a/debian/watch
+++ /dev/null
@@ -1,3 +0,0 @@
-version=4
-
-https://github.com/reedacartwright/ngila/releases .*/archive/v at ANY_VERSION@@ARCHIVE_EXT@
diff --git a/doc/CMakeLists.txt b/doc/CMakeLists.txt
new file mode 100644
index 0000000..e69de29
diff --git a/doc/mkman.pl b/doc/mkman.pl
new file mode 100755
index 0000000..2212ded
--- /dev/null
+++ b/doc/mkman.pl
@@ -0,0 +1,41 @@
+#!/usr/bin/perl -w
+
+use strict;
+
+my $cmds = '';
+while(<>)
+{
+ if(/^\s*XCMD/)
+ {
+ my $m = $_;
+ $m =~ s/^\s*XCMD\(//;
+ my @m = split(/\s*,\s*/, $m);
+ $m[0] =~ s![)][(]!-!g;
+ $m[0] =~ s![()]!!g;
+ $m[1] =~ s![()]!!g;
+ $m[2] =~ s![\"\']!!g;
+ my $arg = ($m[3] =~ /bool/) ? '' : ' arg';
+ my $sm = ($m[1] =~ /\w/) ? " [\\-$m[1]]" : '';
+ $cmds .= ".TP\n.B \\-\\-$m[0]$sm$arg\n$m[2]\n";
+ }
+}
+
+print <<EOF;
+.TH "Ngila" 1
+.SH NAME
+Ngila \- Global Pairwise Alignments with Logarithmic and Affine Gap Costs
+.SH EXAMPLE
+ngila -m zeta -t 0.1 -k 2.0 -r 0.05 -z 1.65 sequences.fasta
+.SH DESCRIPTION
+Ngila is a global, pairwise sequence alignment program that implements the algorithm of
+Miller and Myers (1988), optimized for gap logarithmic and affine gap costs: C(x)=a+b*x+c*ln x. These costs are more biologically realistic and more powerful than standard affine gap costs. In addition Ngila implements two generalized pair-HMM models of indel formation, which transform evolutionary parameters to alignment costs.
+.SH REFERENCE
+Cartwright RA (2007) Ngila: global pairwise alignments with logarithmic and affine gap costs. Bioinformatics. 23(11):1427\-1428
+.SH WEBSITE
+http://scit.us/projects/ngila/
+.SH OPTIONS
+$cmds
+.SH INPUT
+Input sequences must be in FASTA format.
+EOF
+
diff --git a/matrix/CMakeLists.txt b/matrix/CMakeLists.txt
new file mode 100644
index 0000000..e69de29
diff --git a/matrix/dna b/matrix/dna
new file mode 100644
index 0000000..6a5bd68
--- /dev/null
+++ b/matrix/dna
@@ -0,0 +1,5 @@
+ A C G T
+A 0 1 1 1
+C 1 0 1 1
+G 1 1 0 1
+T 1 1 1 0
diff --git a/ngilarc.txt b/ngilarc.txt
new file mode 100644
index 0000000..a008420
--- /dev/null
+++ b/ngilarc.txt
@@ -0,0 +1,4 @@
+model=zeta
+branch-length=0.2
+indel-rate=0.05
+threashold-larger=2000
diff --git a/readme.txt b/readme.txt
new file mode 100644
index 0000000..915e2cc
--- /dev/null
+++ b/readme.txt
@@ -0,0 +1,89 @@
+NGILA VERSION 1.3 - Logarithmic and Affine Sequence Alignments
+
+Copyright (C) 2005-2010 Reed A. Cartwright - All rights reserved.
+
+DESCRIPTION
+ Ngila is a global alignment program that can align pairs of sequences using
+ logarithmic and affine gap penalties.
+
+REFERENCE
+ Cartwright RA (2007) Ngila: global pairwise alignments with logarithmic and
+ affine gap costs. Bioinformatics. 23(11):1427-1428
+
+CONTACT
+ reed at scit.us or racartwright at uh.edu
+
+LICENSE
+ GPL ver 3. See copying.txt.
+
+INSTALLATION
+ Installation from source requires CMake (http://cmake.org/), Boost Libraries
+ (http://boost.org/). Binary packages are available. To install on unix-like
+ sytems simply use
+
+ cmake . && make && make install
+
+ in the extracted source code directory. On Windows you can use CMake GUI
+ to create project files for Visual Studio and install from there. More
+ detailed directions and help can be found on the website.
+
+DOWNLOAD
+ Ngila can be downloaded from the url <http://scit.us/projects/ngila/>, which
+ is its development website.
+
+COMMAND LINE USAGE
+ ngila -m zeta -t 0.1 -k 2.0 -r 0.05 -z 1.65 sequences.fasta
+
+ See 'ngila --help' for complete command line usage.
+
+ See <http://scit.us/projects/ngila/> for more details on running ngila.
+
+MODEL DESCRIPTIONS
+ Ngila includes models of alignment based on evolutionary models. For a basic
+ description of the evolutionary models see Cartwright RA (2009) "Problems and
+ solutions for estimating indel rates and length distributions." Molecular
+ Biology and Evolution, 26:473-480. The models are as follows:
+ zeta: DNA model with indel lengths following a power-law distribution
+ geo: DNA model with a geometric distribution
+ aazeta: protein model (LG 2008) with a power-law distribution
+ aageo: protein model with a geometric distribution
+ cost: specify substitution and gap costs explicitly
+
+INPUT FILES
+ The input file has to be in FASTA or PHYLIP format. If more than two
+ sequences are given then Ngila will align based on the 'pairs' option.
+
+OUTPUT FILES
+ Ngila has two types of output: sequence alignments and distance matrices.
+ Supported sequence alignment formats are Clustal, Fasta, and Phylip. Clustal
+ is the default. The format is read from the output file's extension or
+ specified directly; "ngila -o seqs.fas" and "ngila -o fas:seqs.txt" both
+ produce fasta output. "ngila -o aln:-" sends Clustal formated sequence to
+ stdout. The following extensions are supported for distance matrices:
+ dist-c = likelihood-based cost scores, dist-i = sequence identies, dist-d =
+ sequence distances, dist = lower-triangle like dist-d and upper like dist-i.
+
+SUBSTITUTION MATRIX
+ Used by the "cost" model. An example of the format can be seen in matrix/dna.
+
+NGILARC
+ Commandline options can be specified using an ngilarc file. By default the
+ program looks for $HOME/.ngilarc (unix) or %HOME%/ngilarc.txt (windows). This
+ file can contain long-form command line options like in the ngilarc.txt
+ example file.
+
+
+ALGORITHM
+ Ngila implements a Miller and Myers (1988) candidate list method of sequence
+ alignment with the gap cost being of the form g(x) = a + b*x + c*ln x. Ngila
+ will return the alignment with the minimum cost and has rules for breaking
+ ties. Ngila's main alignment algorithm is divide-and-conquer, which requires
+ O(M) memory; but slower than a holistic, O(MN) memory algorithm.
+
+ Ngila implements a secondary, holistic algorithm for alignment, which is
+ faster. The options -M and -N (-M is for the larger sequence) allow users to
+ specify thresholds for when the holistic algorithm is used instead of the DnC
+ algorithm. For example, command 'ngila -M 5000 -N 5000 seqs.aln' will align
+ the sequences in 'seq.aln' via the divide-and-conquer algorithm, but when
+ subsequences less than or equal to 5000-5000 are being aligned, the holistic
+ algorithm will be used.
diff --git a/releng/build-rel-osx.sh b/releng/build-rel-osx.sh
new file mode 100644
index 0000000..e07b29f
--- /dev/null
+++ b/releng/build-rel-osx.sh
@@ -0,0 +1,55 @@
+#!/bin/sh
+
+PROJ=ngila
+PROJ_DISTS=ngila-1*
+MAKE=make
+CMAKE=cmake
+REPOS=`svn info | grep URL: | perl -pe 's!^URL: (.+)/releng$!$1!'`
+
+echo
+echo Building distributions for $REPOS ...
+
+xcode=/Xcode2.5
+xcodep=${xcode}/usr
+PATH="${xcodep}/bin:${PATH}"; export PATH
+CFLAGS="$CFLAGS -I${xcodep}/include"; export CFLAGS
+CXXFLAGS="$CXXFLAGS -I${xcodep}/include"; export CXXFLAGS
+LDFLAGS="-L${xcodep}/lib"; export LDFLAGS
+CMAKE_OSX_SYSROOT="${xcode}/SDKs/MacOSX10.4u.sdk"; export CMAKE_OSX_SYSROOT
+
+RELENG_DIR=`mktemp -q -d -t ${PROJ}-releng`
+if [ $? -ne 0 ]; then
+ echo "$0: Can't create temp directory, exiting..."
+ exit 1
+fi
+
+echo Using temp directory $RELENG_DIR ...
+echo
+
+DEST_DIR=`pwd`
+SOURCE_DIR="${RELENG_DIR}/source"
+BUILD_DIR="${RELENG_DIR}/build"
+
+svn co -q $REPOS $SOURCE_DIR || exit 1
+
+mkdir $BUILD_DIR || exit 1
+cd $BUILD_DIR || exit 1
+
+$CMAKE $SOURCE_DIR -DCMAKE_BUILD_TYPE=Release \
+ -DCMAKE_OSX_ARCHITECTURES="ppc;i386" \
+ -DCPACK_SYSTEM_NAME=Darwin8-universal \
+ -DAPPLE_BUNDLE=ON
+$MAKE
+$MAKE package
+$MAKE package_source
+
+echo
+echo Moving distribution packages ...
+
+mv -v $PROJ_DISTS $DEST_DIR
+
+echo
+echo Cleaning up ...
+
+cd $DEST_DIR
+rm -rf $RELENG_DIR
diff --git a/releng/build-rel-unix.sh b/releng/build-rel-unix.sh
new file mode 100644
index 0000000..8042aba
--- /dev/null
+++ b/releng/build-rel-unix.sh
@@ -0,0 +1,45 @@
+#!/bin/sh
+
+PROJ=ngila
+PROJ_DISTS=ngila-1*
+MAKE=make
+CMAKE=cmake
+REPOS=`svn info | grep URL: | perl -pe 's!^URL: (.+)/releng$!$1!'`
+
+echo
+echo Building distributions for $REPOS ...
+
+
+RELENG_DIR=`mktemp -q -d -t ${PROJ}-releng.XXX`
+if [ $? -ne 0 ]; then
+ echo "$0: Can't create temp directory, exiting..."
+ exit 1
+fi
+
+echo Using temp directory $RELENG_DIR ...
+echo
+
+DEST_DIR=`pwd`
+SOURCE_DIR="${RELENG_DIR}/source"
+BUILD_DIR="${RELENG_DIR}/build"
+
+svn co -q $REPOS $SOURCE_DIR || exit 1
+
+mkdir $BUILD_DIR || exit 1
+cd $BUILD_DIR || exit 1
+
+$CMAKE $SOURCE_DIR -DCMAKE_BUILD_TYPE=Release -DBoost_USE_STATIC_LIBS=ON -DGSL_USE_STATIC_LIBS=ON
+$MAKE
+$MAKE package
+$MAKE package_source
+
+echo
+echo Moving distribution packages ...
+
+mv -v $PROJ_DISTS $DEST_DIR
+
+echo
+echo Cleaning up ...
+
+cd $DEST_DIR
+rm -rf $RELENG_DIR
diff --git a/releng/build-rel-win.bat b/releng/build-rel-win.bat
new file mode 100644
index 0000000..1171a26
--- /dev/null
+++ b/releng/build-rel-win.bat
@@ -0,0 +1,52 @@
+ at echo off
+
+set PROJ=ngila
+set PROJ_DISTS=ngila-1*
+set MAKE=nmake
+set CMAKE=cmake
+set SVN=svn
+set PERL=perl
+
+%SVN% info | findstr /b URL | %PERL% -pe "s!^URL: (.+)/releng$!$1!" > url.tmp
+set /P REPOS=<url.tmp
+del url.tmp
+
+set PF=%ProgramFiles%
+if defined ProgramFiles(x86) set PF=%ProgramFiles(x86)%
+
+echo.
+echo Building distributions for %REPOS% ...
+
+set RELENG_DIR="%TEMP%\%PROJ%-releng.%RANDOM%"
+mkdir %RELENG_DIR% || exit /B 1
+
+echo Using temp directory %RELENG_DIR% ...
+echo.
+
+set DEST_DIR="%CD%"
+set SOURCE_DIR="%RELENG_DIR%\source"
+set BUILD_DIR="%RELENG_DIR%\build"
+
+%SVN% co -q %REPOS% %SOURCE_DIR% || exit /B 1
+
+mkdir %BUILD_DIR% || exit /B 1
+cd %BUILD_DIR% || exit /B 1
+
+call "%PF%\Microsoft Visual Studio 9.0\VC\vcvarsall.bat" x86
+
+
+%CMAKE% -G "NMake Makefiles" %SOURCE_DIR% -DCMAKE_BUILD_TYPE=Release -DUSE_STATIC_LIBS=ON
+%MAKE%
+%MAKE% package
+%MAKE% package_source
+
+echo.
+echo Copying distribution packages ...
+
+xcopy /Y %PROJ_DISTS% %DEST_DIR%
+
+echo.
+echo Cleaning up ...
+
+cd %DEST_DIR%
+rd /S /Q %RELENG_DIR%
diff --git a/src/CMakeLists.txt b/src/CMakeLists.txt
new file mode 100644
index 0000000..c332248
--- /dev/null
+++ b/src/CMakeLists.txt
@@ -0,0 +1,115 @@
+INCLUDE(CheckIncludeFile)
+INCLUDE(CheckFunctionExists)
+INCLUDE(CheckLibraryExists)
+
+ADD_DEFINITIONS(-DHAVE_CONFIG_H
+ "\"-DPACKAGE_STRING=\\\"${PROJ_NAME} ${PROJ_VERSION}\\\"\""
+ "\"-DPACKAGE_BUGREPORT=\\\"${PROJ_EMAIL}\\\"\""
+ "-DBOOST_ALL_NO_LIB"
+)
+
+SET(NGILA_TABLE_CELL_BITSIZE "16" CACHE STRING
+ "The bitsize of a cell in the travel table. The default allows a maximum gap of 32768 residues.")
+MARK_AS_ADVANCED(NGILA_TABLE_CELL_BITSIZE)
+ADD_DEFINITIONS("-DTABLE_CELL_BITSIZE=${NGILA_TABLE_CELL_BITSIZE}")
+
+IF(USE_THREADS)
+ ADD_DEFINITIONS("-DUSE_THREADS=1")
+ SET(USING_BOOST_THREAD "thread")
+ENDIF(USE_THREADS)
+
+IF(CMAKE_COMPILER_IS_GNUCC)
+ SET(CMAKE_C_FLAGS "${CMAKE_C_FLAGS} -Wall")
+ENDIF(CMAKE_COMPILER_IS_GNUCC)
+
+IF(CMAKE_COMPILER_IS_GNUCC)
+ SET(CMAKE_C_FLAGS "${CMAKE_C_FLAGS} -Wall")
+ENDIF(CMAKE_COMPILER_IS_GNUCC)
+
+IF(MSVC)
+ ADD_DEFINITIONS(-D_CRT_SECURE_NO_DEPRECATE -D_CRT_NONSTDC_NO_DEPRECATE)
+ ADD_DEFINITIONS(-D_SCL_SECURE_NO_DEPRECATE)
+ SET(CMAKE_C_FLAGS "${MY_C_FLAGS} /wd4100 /wd4101 /wd4189 /wd4706")
+ IF(USE_STATIC_LIBS)
+ STRING(REGEX REPLACE "/MD" "/MT" CMAKE_C_FLAGS_RELEASE "${CMAKE_C_FLAGS_RELEASE}")
+ STRING(REGEX REPLACE "/MDd" "/MTd" CMAKE_C_FLAGS_DEBUG "${CMAKE_C_FLAGS_DEBUG}")
+ STRING(REGEX REPLACE "/MD" "/MT" CMAKE_CXX_FLAGS_RELEASE "${CMAKE_CXX_FLAGS_RELEASE}")
+ STRING(REGEX REPLACE "/MDd" "/MTd" CMAKE_CXX_FLAGS_DEBUG "${CMAKE_CXX_FLAGS_DEBUG}")
+ ENDIF(USE_STATIC_LIBS)
+
+ # Use the highest warning level for visual studio.
+ SET(CMAKE_C_WARNING_LEVEL 4)
+ IF(CMAKE_C_FLAGS MATCHES "/W[0-4]")
+ STRING(REGEX REPLACE "/W[0-4]" "/W4"
+ CMAKE_C_FLAGS "${CMAKE_C_FLAGS}")
+ ELSE(CMAKE_C_FLAGS MATCHES "/W[0-4]")
+ SET(CMAKE_C_FLAGS "${CMAKE_C_FLAGS} /W4")
+ ENDIF(CMAKE_C_FLAGS MATCHES "/W[0-4]")
+ENDIF(MSVC)
+
+IF(WIN32)
+ SET(RES_SRCS ${RES_SRCS} ngila.rc)
+ ADD_DEFINITIONS(-DNOMINMAX)
+ENDIF(WIN32)
+
+IF(APPLE_BUNDLE_APP)
+ SET(PROJ_MACOSX_BUNDLE MACOSX_BUNDLE)
+ENDIF(APPLE_BUNDLE_APP)
+
+SET(Boost_ADDITIONAL_VERSIONS 1.50 1.49 1.48 1.47 1.46 1.45 1.44 1.43 1.42 1.41 1.40)
+
+FIND_PACKAGE(Boost 1.36.0 REQUIRED COMPONENTS program_options ${USING_BOOST_THREAD})
+IF(Boost_FOUND)
+ ADD_DEFINITIONS(-DBOOST_ALL_NO_LIB -DBOOST_PROGRAM_OPTIONS_NO_LIB -DBOOST_TYPEOF_SILENT)
+ INCLUDE_DIRECTORIES(${Boost_INCLUDE_DIRS})
+ LINK_DIRECTORIES(${Boost_LIBRARY_DIRS})
+ IF(NOT Boost_USE_STATIC_LIBS)
+ ADD_DEFINITIONS(-DBOOST_DYN_LINK -DBOOST_PROGRAM_OPTIONS_DYN_LINK)
+ ENDIF(NOT Boost_USE_STATIC_LIBS)
+ #SET(CMAKE_CXX_FLAGS_DEBUG "${CMAKE_CXX_FLAGS_DEBUG} ${Boost_LIB_DIAGNOSTIC_DEFINITIONS}")
+ENDIF(Boost_FOUND)
+
+FIND_PACKAGE(GSL)
+IF(GSL_FOUND)
+ INCLUDE_DIRECTORIES(${GSL_INCLUDE_DIR})
+ LINK_DIRECTORIES(${GSL_LINK_DIRECTORIES})
+ SET(CMAKE_EXE_LINKER_FLAGS "${CMAKE_EXE_LINKER_FLAGS} ${GSL_EXE_LINKER_FLAGS}")
+ ADD_DEFINITIONS(-DUSE_GSL_ZETA)
+ENDIF(GSL_FOUND)
+
+CONFIGURE_FILE(${CMAKE_CURRENT_SOURCE_DIR}/config.h.cmake ${CMAKE_CURRENT_BINARY_DIR}/config.h)
+INCLUDE_DIRECTORIES(BEFORE "${CMAKE_CURRENT_BINARY_DIR}")
+
+ADD_EXECUTABLE(ngila ${PROJ_MACOSX_BUNDLE}
+ align.cpp matparser.cpp models.cpp ngila.cpp seqdb.cpp
+ align.h matparser.h models.h ngila.h ngila_app.h seqdb.h seqparser.h
+ xm.h ngila.cmds sort.h
+ "${CMAKE_CURRENT_BINARY_DIR}/config.h"
+ ${COMPAT}
+ ${RES_SRCS})
+
+IF(USE_STATIC_LIBS)
+ SET_TARGET_PROPERTIES(ngila PROPERTIES LINK_SEARCH_END_STATIC On)
+ENDIF(USE_STATIC_LIBS)
+
+TARGET_LINK_LIBRARIES(ngila ${CMAKE_REQUIRED_LIBRARIES})
+TARGET_LINK_LIBRARIES(ngila ${Boost_LIBRARIES})
+TARGET_LINK_LIBRARIES(ngila ${GSL_LIBRARIES})
+
+IF(APPLE_BUNDLE_APP)
+ SET_TARGET_PROPERTIES(ngila PROPERTIES OUTPUT_NAME ${PROJ_BUNDLE_NAME})
+ SET(PROJ_INSTALL_EXTRA_ARGS BUNDLE DESTINATION "${PROJ_BUNDLE_LOCATION}")
+ENDIF(APPLE_BUNDLE_APP)
+
+INSTALL(TARGETS ngila ${PROJ_INSTALL_EXTRA_ARGS} RUNTIME DESTINATION bin)
+
+IF(UNIX)
+ INSTALL(FILES ngila.desktop DESTINATION share/applications )
+ INSTALL(FILES ngila.png DESTINATION share/pixmaps )
+ENDIF(UNIX)
+
+IF(APPLE_BUNDLE_APP)
+ INSTALL(CODE "file(MAKE_DIRECTORY \$ENV{DESTDIR}\${CMAKE_INSTALL_PREFIX}/bin)")
+ INSTALL(CODE "execute_process(COMMAND ln -s \"../MacOS/${PROJ_BUNDLE_NAME}\" ngila
+ WORKING_DIRECTORY \$ENV{DESTDIR}\${CMAKE_INSTALL_PREFIX}/bin)")
+ENDIF(APPLE_BUNDLE_APP)
diff --git a/src/align.cpp b/src/align.cpp
new file mode 100644
index 0000000..47ee9d0
--- /dev/null
+++ b/src/align.cpp
@@ -0,0 +1,596 @@
+/****************************************************************************
+ * Copyright (C) 2005-2010 Reed A. Cartwright, PhD <reed at scit.us> *
+ * *
+ * This program is free software: you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation, either version 3 of the License, or *
+ * (at your option) any later version. *
+ * *
+ * This program is distributed in the hope that it will be useful, *
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of *
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
+ * GNU General Public License for more details. *
+ * *
+ * You should have received a copy of the GNU General Public License *
+ * along with this program. If not, see <http://www.gnu.org/licenses/>. *
+ ****************************************************************************/
+
+#ifdef _MSC_VER
+# pragma warning(disable: 4288)
+#endif
+
+#include <algorithm>
+
+#include "ngila.h"
+#include "align.h"
+
+using namespace std;
+
+template<class T>
+inline T min3(T a, T b, T c)
+{
+ return (std::min)(a, (std::min)(b,c));
+}
+
+inline bool dle(double a, double b)
+{
+ return (a <= b*(1.0+DBL_EPSILON));
+}
+
+double aligner::align(alignment &aln)
+{
+ aln.data.clear();
+ return align_x(aln.seqA.dna, aln.seqB.dna, aln.data);
+}
+
+double aligner::align_x(const sequence &seqA, const sequence &seqB, aln_data &rAln)
+{
+ size_t sz = seqB.size()+1;
+ CC[0].resize(sz);
+ CC[1].resize(sz);
+ RR[0].resize(sz);
+ RR[1].resize(sz);
+ SF.resize(sz);
+ SR.resize(sz);
+ DM.resize(sz);
+ GC.resize(1u+max(seqA.size(),seqB.size()));
+ FGC.resize(GC.size());
+ for(size_t u = 0;u<GC.size();++u)
+ {
+ GC[u] = costs.gapcost(u);
+ FGC[u] = costs.freegapcost(u);
+ }
+
+ // Allocate travel table to maximum possible size
+ size_t szA = (std::min)(szMa,seqA.size())+1;
+ tabTravel.resize(szA);
+ size_t szB = (std::min)(szMb,seqB.size())+1;
+ for(travel_table::iterator it = tabTravel.begin();
+ it != tabTravel.end(); ++it)
+ it->resize(szB);
+
+ alnBuf.reserve(szA+szB);
+
+ //run recursive algorithm
+ return align_r(seqA.begin(), seqA.end(), seqB.begin(), seqB.end(), rAln, bFreeEnds, bFreeEnds);
+}
+
+double aligner::align_mn(sequence::const_iterator itA1, sequence::const_iterator itA2,
+ sequence::const_iterator itB1, sequence::const_iterator itB2,
+ aln_data &rAln, bool bFreeFront, bool bFreeBack)
+{
+ // O(MN) memory algorithm
+ size_t szNa = itA2-itA1;
+ size_t szNb = itB2-itB1;
+
+ CC[0][0] = 0.0;
+ tabTravel[0][0] = 0;
+ for(size_t j=1;j<=szNb;++j)
+ {
+ CC[0][j] = GC[j];
+ tabTravel[0][j] = static_cast<travel_cell>(j);
+ }
+ // TODO: Is this needed?
+ SF[0].assign(1, indel(0, szNa, CC[0][0]));
+
+ for(size_t i=1;i<=szNa;++i)
+ {
+ CC[1][0] = (bFreeFront ? FGC : GC)[i];
+ tabTravel[i][0] = -static_cast<travel_cell>(i);
+ for(size_t j=1;j<=szNb;++j)
+ {
+ double dM = CC[0][j-1]+costs.mCost[(size_t)itA1[i-1]][(size_t)itB1[j-1]];
+ update_ins_forward(T,i,j,szNb);
+ double dI = indel_cost(T.back(),j);
+ double dD;
+ if(bFreeBack && j == szNb)
+ {
+ update_del_forward_f(SF[j],i,j,szNa);
+ dD = indel_fcost(SF[j].back(),i);
+ }
+ else
+ {
+ update_del_forward(SF[j],i,j,szNa);
+ dD = indel_cost(SF[j].back(),i);
+ }
+
+ if(dM < dI && dM < dD)
+ {
+ CC[1][j] = dM;
+ tabTravel[i][j] = static_cast<travel_cell>(0);
+ }
+ else if(dI <= dM && dI < dD)
+ {
+ CC[1][j] = dI;
+ tabTravel[i][j] = static_cast<travel_cell>(j-T.back().p);
+ }
+ else
+ {
+ CC[1][j] = dD;
+ tabTravel[i][j] = -static_cast<travel_cell>(i-SF[j].back().p);
+ }
+ }
+ swap(CC[0], CC[1]);
+ }
+ size_t i = szNa, j = szNb;
+
+ alnBuf.clear();
+ while(!(i == 0 && j == 0))
+ {
+ travel_cell t = tabTravel[i][j];
+ alnBuf.push_back(t);
+ if(t == 0)
+ {
+ --i;
+ --j;
+ }
+ else if(t > 0)
+ j -= static_cast<size_t>(t);
+ else
+ i -= static_cast<size_t>(-t);
+ }
+ rAln.insert(rAln.end(), alnBuf.rbegin(), alnBuf.rend());
+ return CC[0][szNb];
+}
+
+double aligner::align_s(sequence::const_iterator itA1, sequence::const_iterator itA2,
+ sequence::const_iterator itB1, sequence::const_iterator itB2,
+ aln_data &rAln, bool bFreeFront, bool bFreeBack)
+{
+ size_t szNa = itA2-itA1;
+ size_t szNb = itB2-itB1;
+ if(szNa == 0 && szNb == 0)
+ return 0.0; // Nothing to do
+ else if(szNb == 0)
+ {
+ rAln.push_back(-static_cast<alignment::aln_atom>(szNa));
+ return (bFreeFront||bFreeBack ? FGC : GC)[szNa]; // delete A
+ }
+ else if(szNa == 0)
+ {
+ rAln.push_back(static_cast<alignment::aln_atom>(szNb));
+ return GC[szNb]; // insert B
+ }
+ else if(szNa == 1 && szNb == 1)
+ {
+ double d1 = costs.mCost[(size_t)itA1[0]][(size_t)itB1[0]];
+ double d2 = (bFreeFront||bFreeBack) ? FGC[1]+GC[1] : GC[1]+GC[1];
+ if(d1 <= d2)
+ {
+ rAln.push_back(0);
+ return d1;
+ }
+ else if(bFreeFront)
+ {
+ rAln.push_back(-1);
+ rAln.push_back(1);
+ }
+ else
+ {
+ rAln.push_back(1);
+ rAln.push_back(-1);
+ }
+ return d2;
+ }
+ else if(szNa == 1)
+ {
+ double dTemp = costs.mCost[(size_t)itA1[0]][(size_t)itB1[0]] + GC[szNb-1];
+ double dMin = dTemp;
+ size_t i = 0;
+ for(size_t j=1;j<szNb-1;++j)
+ {
+ dTemp = GC[j] + costs.mCost[(size_t)itA1[0]][(size_t)itB1[j]] + GC[szNb-j-1];
+ if(dTemp < dMin)
+ {
+ dMin = dTemp;
+ i = j;
+ }
+ }
+ dTemp = costs.mCost[(size_t)itA1[0]][(size_t)itB1[szNb-1]] + GC[szNb-1];
+ if(dTemp < dMin)
+ {
+ dMin = dTemp;
+ i = szNb-1;
+ }
+ dTemp = (bFreeFront ? FGC : GC)[1] + GC[szNb];
+ if(dTemp < dTemp)
+ {
+ dMin = dTemp;
+ i = static_cast<size_t>(-1);
+ }
+
+ dTemp = (bFreeBack ? FGC : GC)[1] + GC[szNb];
+ if(dTemp < dMin)
+ {
+ dMin = dTemp;
+ i = szNb;
+ }
+
+ if(i == static_cast<size_t>(-1))
+ {
+ // Del A, Ins B(1,Nb)
+ rAln.push_back(-1);
+ rAln.push_back(static_cast<alignment::aln_atom>(szNb));
+ }
+ else if( i == 0)
+ {
+ // A->B(1), Ins B(2,Nb)
+ rAln.push_back(0);
+ rAln.push_back(static_cast<alignment::aln_atom>(szNb-1));
+ }
+ else if(i == szNb-1)
+ {
+ // Ins B(1,N-1), A->B(Nb)
+ rAln.push_back(static_cast<alignment::aln_atom>(szNb-1));
+ rAln.push_back(0);
+ }
+ else if( i == szNb)
+ {
+ // Ins B(1,N), Del A
+ rAln.push_back(static_cast<alignment::aln_atom>(szNb));
+ rAln.push_back(-1);
+ }
+ else
+ {
+ //Ins B(1,i), A->B(i+1), Ins B(i+2,Nb)
+ rAln.push_back(static_cast<alignment::aln_atom>(i));
+ rAln.push_back(0);
+ rAln.push_back(static_cast<alignment::aln_atom>(szNb-i-1));
+ }
+ return dMin;
+ }
+ else if(szNb == 1)
+ {
+
+ double dTemp = costs.mCost[(size_t)itA1[0]][(size_t)itB1[0]]
+ + (bFreeBack ? FGC : GC)[szNa-1];
+ double dMin = dTemp;
+ size_t i = 0;
+ for(size_t j=1;j<szNa-1;++j)
+ {
+ dTemp = costs.mCost[(size_t)itA1[j]][(size_t)itB1[0]]
+ + (bFreeFront ? FGC : GC)[j]
+ + (bFreeBack ? FGC : GC)[szNa-j-1];
+ if(dTemp < dMin)
+ {
+ dMin = dTemp;
+ i = j;
+ }
+ }
+ dTemp = costs.mCost[(size_t)itA1[szNa-1]][(size_t)itB1[0]]
+ + (bFreeFront ? FGC : GC)[szNa-1];
+ if(dTemp < dMin)
+ {
+ dMin = dTemp;
+ i = szNa-1;
+ }
+
+ dTemp = GC[1] + (bFreeBack ? FGC : GC)[szNa];
+ if(dTemp < dMin)
+ {
+ dMin = dTemp;
+ i = (size_t)-1;
+ }
+ dTemp = (bFreeFront ? FGC : GC)[szNa] + GC[1];
+ if(dTemp < dMin)
+ {
+ dMin = dTemp;
+ i = szNa;
+ }
+
+ if(i == (size_t)-1)
+ {
+ // Ins B, Del A(1,Na)
+ rAln.push_back(1);
+ rAln.push_back(-static_cast<alignment::aln_atom>(szNa));
+ }
+ else if(i == 0)
+ {
+ // B->A(1), Del A(2,Na)
+ rAln.push_back(0);
+ rAln.push_back(-static_cast<alignment::aln_atom>(szNa-1));
+ }
+ else if(i == szNa-1)
+ {
+ // Del A(1,Na-1), B->A(Na)
+ rAln.push_back(-static_cast<alignment::aln_atom>(szNa-1));
+ rAln.push_back(0);
+ }
+ else if(i == szNa)
+ {
+ // Del A(1,Na), Ins B
+ rAln.push_back(-static_cast<alignment::aln_atom>(szNa));
+ rAln.push_back(1);
+ }
+ else
+ {
+ //Del A(1,i), B->A(i+1), Del A(i+2,Na)
+ rAln.push_back(-static_cast<alignment::aln_atom>(i));
+ rAln.push_back(0);
+ rAln.push_back(-static_cast<alignment::aln_atom>(szNa-i-1));
+ }
+ return dMin;
+ }
+ return 0.0;
+}
+
+
+double aligner::align_r(sequence::const_iterator itA1, sequence::const_iterator itA2,
+ sequence::const_iterator itB1, sequence::const_iterator itB2,
+ aln_data &rAln, bool bFreeFront, bool bFreeBack)
+{
+ size_t szNa = itA2-itA1;
+ size_t szNb = itB2-itB1;
+ size_t szMh = szNa/2;
+
+ if(szNa <= 1 || szNb <= 1) // Handle Special Cases
+ return align_s(itA1, itA2, itB1, itB2, rAln, bFreeFront, bFreeBack);
+ else if(szNa <= szMa && szNb <= szMb)
+ return align_mn(itA1, itA2, itB1, itB2, rAln, bFreeFront, bFreeBack);
+
+ CC[0][0] = 0.0;
+ for(size_t j=1;j<=szNb;++j)
+ CC[0][j] = GC[j];
+
+ SF[0].assign(1, indel(0, szNa, CC[0][0]));
+ // Forward Algorithm
+ for(size_t i=1;i<=szMh;++i)
+ {
+ CC[1][0] = (bFreeFront ? FGC : GC)[i];
+ for(size_t j=1;j<=szNb;++j)
+ {
+ update_ins_forward(T,i,j,szNb);
+ update_del_forward(SF[j],i,j,szNa);
+ CC[1][j] = min3(CC[0][j-1]+costs.mCost[(size_t)itA1[i-1]][(size_t)itB1[j-1]],
+ indel_cost(T.back(),j), indel_cost(SF[j].back(),i));
+ }
+ swap(CC[0], CC[1]);
+ }
+ //CC[0] is now C(M/2,j)
+
+ //Reverse Algorithm
+ RR[0][szNb] = 0.0;
+ DM[szNb].c = bFreeBack ? indel_fcost(SF[szNb][0],szNa)
+ : indel_cost(SF[szNb][0],szNa);
+ DM[szNb].s = 0;
+ DM[szNb].z = 0;
+ DM[szNb].x = szNa;
+ RR[0][0] = GC[szNb];
+ DM[0].c = indel_cost(SF[0][0],szNa)+RR[0][0];
+ DM[0].s = 0;
+ DM[0].z = 0;
+ DM[0].x = szNa;
+ for(size_t j=szNb-1;j>0;--j)
+ {
+ RR[0][j] = GC[szNb-j];
+ DM[j].c = indel_cost(SF[j][0],szNa)+RR[0][j];
+ DM[j].s = 0;
+ DM[j].z = 0;
+ DM[j].x = szNa;
+ }
+ for(size_t i=szNa-1;i!=szMh-1;--i)
+ {
+ if( SF[szNb].size() < DM[szNb].z+1 && i <= SF[szNb][DM[szNb].z+1].x)
+ ++DM[szNb].z; // Advance position
+ RR[1][szNb] = (bFreeBack ? FGC : GC)[szNa-i];
+ double dTemp = indel_cost(SF[szNb][DM[szNb].z],i)+RR[1][szNb];
+ if(dTemp < DM[szNb].c)
+ {
+ DM[szNb].c = dTemp;
+ DM[szNb].s = DM[szNb].z;
+ DM[szNb].x = i;
+ }
+
+ for(size_t j=szNb-1;j!=(size_t)-1;--j)
+ {
+ update_ins_reverse(T,i,j,szNb);
+ update_del_reverse(SR[j],i,j,szNa);
+ // Type I
+ double d1 = RR[0][j+1]+costs.mCost[(size_t)itA1[i]][(size_t)itB1[j]];
+ double d2 = indel_rcost(T.back(),j);
+ double d3 = indel_rcost(SR[j].back(),i);
+ RR[1][j] = min3(d1,d2,d3);
+
+ // Minimum Type II cost
+ if( SF[j].size() < DM[j].z+1 && i <= SF[j][DM[j].z+1].x)
+ ++DM[j].z; // Advance position
+ if(bFreeFront && j==0)
+ dTemp = indel_fcost(SF[j][DM[j].z],i)+RR[1][j];
+ else
+ dTemp = indel_cost(SF[j][DM[j].z],i)+RR[1][j];
+ if(dTemp < DM[j].c)
+ {
+ DM[j].c = dTemp;
+ DM[j].s = DM[j].z;
+ DM[j].x = i;
+ }
+ }
+ swap(RR[0], RR[1]);
+ }
+
+ double dMin = CC[0][0]+RR[0][0];
+ size_t pp = szMh;
+ size_t xx = pp;
+ size_t jj = 0;
+ if(dle(DM[0].c, dMin))
+ {
+ dMin = DM[0].c;
+ pp = SF[0][DM[0].s].p;
+ xx = DM[0].x;
+ }
+ for(size_t j=1;j<=szNb;++j)
+ {
+ double dTemp = CC[0][j]+RR[0][j];
+ if(dTemp < dMin)
+ {
+ dMin = dTemp;
+ pp = szMh;
+ xx = pp;
+ jj = j;
+ }
+ dTemp = DM[j].c;
+ if(dle(dTemp,dMin))
+ {
+ dMin = dTemp;
+ pp = SF[j][DM[j].s].p;
+ xx = DM[j].x;
+ jj = j;
+ }
+ }
+
+ align_r(itA1, itA1+pp, itB1, itB1+jj, rAln, bFreeFront, false);
+ if(xx != pp)
+ rAln.push_back(-static_cast<alignment::aln_atom>(xx-pp)); // Delete itA1+pp .. itA1+xx
+ align_r(itA1+xx, itA2, itB1+jj, itB2, rAln, false, bFreeBack);
+ return dMin;
+}
+
+void aligner::update_ins_forward(indel_vec &T, size_t /*i*/, size_t j, size_t szZ)
+{
+ if(j == 1)
+ {
+ T.clear();
+ T.push_back(indel(0, szZ, CC[1][0]));
+ return;
+ }
+ if( j > T.back().x)
+ T.pop_back();
+ if( CC[1][j-1] + GC[1] < indel_cost(T.back(), j) ||
+ CC[1][j-1] + GC[szZ-j+1] <= indel_cost(T.back(), szZ) )
+ {
+ while(T.size() && (CC[1][j-1] + GC[T.back().x - j+1] < indel_cost_x(T.back()) ||
+ CC[1][j-1] + GC[szZ-j+1] <= indel_cost(T.back(), szZ)))
+ T.pop_back();
+ T.push_back(indel(j-1, (T.size() ?
+ (T.back().p+costs.kstar(j-1-T.back().p, CC[1][j-1]-T.back().d))
+ : szZ) ,CC[1][j-1]));
+ }
+}
+
+void aligner::update_del_forward(indel_vec &T, size_t i, size_t j, size_t szZ)
+{
+ if(i == 1)
+ {
+ T.clear();
+ T.push_back(indel(0, szZ, CC[0][j]));
+ return;
+ }
+ if( i > T.back().x)
+ T.pop_back();
+ if( CC[0][j] + GC[1] < indel_cost(T.back(), i) ||
+ CC[0][j] + GC[szZ-i+1] <= indel_cost(T.back(),szZ) )
+ {
+ while(T.size() && (CC[0][j] + GC[T.back().x - i+1] < indel_cost_x(T.back()) ||
+ CC[0][j] + GC[szZ-i+1] <= indel_cost(T.back(), szZ)))
+ T.pop_back();
+ T.push_back(indel(i-1, (T.size() ?
+ (T.back().p+costs.kstar(i-1-T.back().p, CC[0][j]-T.back().d))
+ : szZ), CC[0][j]));
+ }
+}
+
+void aligner::update_ins_reverse(indel_vec &T, size_t i, size_t j, size_t szZ)
+{
+ if(j == szZ-1)
+ {
+ T.clear();
+ T.push_back(indel(szZ, 0, RR[1][szZ]));
+ return;
+ }
+ if( j < T.back().x)
+ T.pop_back();
+ if( RR[1][j+1] + GC[1] < indel_rcost(T.back(),j) ||
+ RR[1][j+1] + GC[j+1] <= indel_rcost(T.back(),0) )
+ {
+ while(T.size() && (RR[1][j+1] + GC[j+1-T.back().x] < indel_rcost_x(T.back()) ||
+ RR[1][j+1] + GC[j+1] <= indel_rcost(T.back(),0)))
+ T.pop_back();
+ T.push_back(indel(j+1, (T.size() ?
+ (T.back().p-costs.kstar(T.back().p-(j+1), RR[1][j+1]-T.back().d))
+ : 0) ,RR[1][j+1]));
+ }
+}
+
+void aligner::update_del_forward_f(indel_vec &T, size_t i, size_t j, size_t szZ)
+{
+ if(i == 1)
+ {
+ T.clear();
+ T.push_back(indel(0, szZ, CC[0][j]));
+ return;
+ }
+ if( i > T.back().x)
+ T.pop_back();
+ if( CC[0][j] + FGC[1] < indel_fcost(T.back(), i) ||
+ CC[0][j] + FGC[szZ-i+1] <= indel_fcost(T.back(),szZ) )
+ {
+ while(T.size() && (CC[0][j] + FGC[T.back().x - i+1] < indel_fcost_x(T.back()) ||
+ CC[0][j] + FGC[szZ-i+1] <= indel_fcost(T.back(), szZ)))
+ T.pop_back();
+ T.push_back(indel(i-1, (T.size() ?
+ (T.back().p+costs.kstar_f(i-1-T.back().p, CC[0][j]-T.back().d))
+ : szZ), CC[0][j]));
+ }
+}
+
+void aligner::update_del_reverse(indel_vec &T, size_t i, size_t j, size_t szZ)
+{
+ if(i == szZ-1)
+ {
+ T.clear();
+ T.push_back(indel(szZ, 0, RR[0][j]));
+ return;
+ }
+ if( i < T.back().x)
+ T.pop_back();
+ if( RR[0][j] + GC[1] < indel_rcost(T.back(),i) ||
+ RR[0][j] + GC[i+1] <= indel_rcost(T.back(),0) )
+ {
+ while(T.size() && (RR[0][j] + GC[i+1-T.back().x] < indel_rcost_x(T.back()) ||
+ RR[0][j] + GC[i+1] <= indel_rcost(T.back(),0)))
+ T.pop_back();
+ T.push_back(indel(i+1, (T.size() ?
+ (T.back().p-costs.kstar(T.back().p-(i+1), RR[0][j]-T.back().d))
+ : 0) ,RR[0][j]));
+ }
+}
+
+double alignment::identity() const {
+ aln_data::size_type matches = 0, mismatches = 0;
+ std::string::const_iterator ait = seqA.dna.begin(), bit = seqB.dna.begin();
+ for(aln_data::const_iterator cit = data.begin(); cit != data.end(); ++cit) {
+ if(*cit == 0) {
+ if(*ait == *bit)
+ ++matches;
+ else
+ ++mismatches;
+ ++ait;
+ ++bit;
+ } else if(*cit > 0)
+ bit+=*cit;
+ else
+ ait-=*cit;
+ }
+ return (1.0*matches)/(matches+mismatches);
+}
+
diff --git a/src/align.h b/src/align.h
new file mode 100644
index 0000000..d9a9ec7
--- /dev/null
+++ b/src/align.h
@@ -0,0 +1,235 @@
+/****************************************************************************
+ * Copyright (C) 2005-2010 Reed A. Cartwright, PhD <reed at scit.us> *
+ * *
+ * This program is free software: you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation, either version 3 of the License, or *
+ * (at your option) any later version. *
+ * *
+ * This program is distributed in the hope that it will be useful, *
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of *
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
+ * GNU General Public License for more details. *
+ * *
+ * You should have received a copy of the GNU General Public License *
+ * along with this program. If not, see <http://www.gnu.org/licenses/>. *
+ ****************************************************************************/
+
+#ifndef ALIGN_H
+#define ALIGN_H
+
+#include <vector>
+#include <iostream>
+#include <iomanip>
+#include <ios>
+
+#include "seqdb.h"
+#include "models.h"
+
+const char chGap = '-';
+
+class aligner;
+
+class alignment
+{
+public:
+ typedef int aln_atom;
+ typedef std::vector<aln_atom> aln_data;
+
+ alignment(const seq_data &A, const seq_data &B) : seqA(A), seqB(B) { }
+
+ template<class OS>
+ inline void print(OS &os, int format=0, double cost=0.0, int dir=0, bool swapped=false) const;
+
+ double identity() const;
+
+protected:
+ const seq_data &seqA, &seqB;
+
+ aln_data data;
+
+ friend class aligner;
+};
+
+class aligner {
+public:
+ struct indel {
+ indel() { }
+ indel(size_t pp, size_t xx, double dd)
+ : p(pp), x(xx), d(dd) { }
+ size_t p; // Column or Row
+ size_t x; // Crossing Point
+ double d; // Score
+ };
+ struct min_mid_cost {
+ double c; // Minimum cost
+ size_t s; // minimum position in s-vector
+ size_t x; // bottom of minimum indel
+ size_t z; // active position in s-vector
+ };
+
+ typedef std::string sequence;
+ typedef alignment::aln_data aln_data;
+ typedef std::vector<indel> indel_vec;
+
+
+ typedef _travel_cell travel_cell;
+ typedef std::vector<travel_cell> travel_row;
+ typedef std::vector<travel_row> travel_table;
+
+ aligner(const cost_model& m, size_t ma, size_t mb, bool fe) : costs(m), szMa(ma), szMb(mb), bFreeEnds(fe) { };
+
+ double align(alignment &rAln);
+
+protected:
+ cost_model costs;
+ size_t szMa;
+ size_t szMb;
+ bool bFreeEnds;
+
+ travel_table tabTravel;
+ aln_data alnBuf;
+
+private:
+ inline double indel_cost(const indel& in, size_t q) const { return in.d + GC[q-in.p]; }
+ inline double indel_cost_x(const indel& in) const { return in.d + GC[in.x-in.p]; }
+ inline double indel_fcost(const indel& in, size_t q) const { return in.d + FGC[q-in.p]; }
+ inline double indel_fcost_x(const indel& in) const { return in.d + FGC[in.x-in.p]; }
+ inline double indel_rcost(const indel& in, size_t q) const { return in.d + GC[in.p-q]; }
+ inline double indel_rcost_x(const indel& in) const { return in.d + GC[in.p-in.x]; }
+
+ double align_x(const sequence &seqa, const sequence &seqB, aln_data &rAln);
+
+ double align_mn(sequence::const_iterator itA1, sequence::const_iterator itA2,
+ sequence::const_iterator itB1, sequence::const_iterator itB2,
+ aln_data &rAln, bool bFreeFront, bool bFreeBack);
+
+ double align_r(sequence::const_iterator itA1, sequence::const_iterator itA2,
+ sequence::const_iterator itB1, sequence::const_iterator itB2,
+ aln_data &rAln, bool bFreeFront, bool bFreeBack);
+
+ double align_s(sequence::const_iterator itA1, sequence::const_iterator itA2,
+ sequence::const_iterator itB1, sequence::const_iterator itB2,
+ aln_data &rAln, bool bFreeFront, bool bFreeBack);
+
+
+ void update_ins_forward(indel_vec &T, size_t i, size_t j, size_t szZ);
+ void update_del_forward(indel_vec &T, size_t i, size_t j, size_t szZ);
+ void update_del_forward_f(indel_vec &T, size_t i, size_t j, size_t szZ);
+ void update_ins_reverse(indel_vec &T, size_t i, size_t j, size_t szZ);
+ void update_del_reverse(indel_vec &T, size_t i, size_t j, size_t szZ);
+
+ std::vector<double> CC[2];
+ std::vector<double> RR[2];
+ std::vector<double> GC, FGC;
+
+ indel_vec T;
+ std::vector<indel_vec> SF;
+ std::vector<indel_vec> SR;
+
+ std::vector<min_mid_cost> DM;
+};
+
+template<class OS>
+inline void alignment::print(OS &os, int format, double cost, int dir, bool swapped) const
+{
+ // 0 = clustal, 1 = fasta, 2 = phylip
+ std::string strA, strB, strC, nameA, nameB;
+ if(format == 0) {
+ nameA = seqA.name.substr(0, 14);
+ nameB = seqB.name.substr(0, 14);
+ } else if(format == 2) {
+ nameA = seqA.name.substr(0, 9);
+ nameB = seqB.name.substr(0, 9);
+ } else {
+ nameA = seqA.name;
+ nameB = seqB.name;
+ }
+ //preallocate
+ strA.reserve(4048);
+ strB.reserve(4048);
+ strC.reserve(4048);
+
+ std::string::const_iterator ait = seqA.dna.begin(), bit = seqB.dna.begin();
+ for(aln_data::const_iterator cit = data.begin(); cit != data.end(); ++cit) {
+ if(*cit == 0) {
+ if(*ait == *bit)
+ strC.append(1, '*');
+ else
+ strC.append(1, ' ');
+ strA.append(ait, ait+1); ++ait;
+ strB.append(bit, bit+1); ++bit;
+ } else if(*cit > 0) {
+ strC.append(*cit, ' ');
+ strA.append(*cit, chGap);
+ strB.append(bit, bit+*cit); bit+=*cit;
+ } else {
+ strC.append(-*cit, ' ');
+ strA.append(ait, ait-*cit); ait-=*cit;
+ strB.append(-*cit, chGap);
+ }
+ }
+ if(swapped) {
+ swap(strA, strB);
+ swap(nameA, nameB);
+ }
+ seq_transform(strA, dir);
+ seq_transform(strB, dir);
+ seq_transform(strC, dir);
+
+ size_t sz = strA.size();
+ std::string ss;
+ if(format == 0) {
+ os << "CLUSTAL multiple sequence alignment (Created by " << PACKAGE_STRING
+ << "; Cost = " << std::setprecision(10) << cost
+ << ")\n\n";
+
+ // Print interleaved sequences
+ size_t a=0, b=0;
+ for(size_t u = 0; u < sz; u+=60) {
+ // Print a row of each sequence
+ ss = strA.substr(u, 60);
+ a += ss.length() - std::count(ss.begin(), ss.end(), chGap);
+ os << std::setw(14) << std::setiosflags(std::ios::left) << nameA
+ << " " << ss << " " << a << "\n";
+
+ ss = strB.substr(u, 60);
+ b += ss.length() - std::count(ss.begin(), ss.end(), chGap);
+ os << std::setw(14) << std::setiosflags(std::ios::left) << nameB
+ << " " << ss << " " << b << "\n";
+
+ os << std::setw(15) << std::setiosflags(std::ios::left)
+ << " " << strC.substr(u, 60) << "\n\n";
+ }
+ } else if(format == 1) {
+ os << ">" << nameA << " "
+ << "(Created by " << PACKAGE_STRING
+ << "; Cost = " << std::setprecision(10) << cost
+ << ")\n";
+ for(size_t u = 0; u < sz; u += 80)
+ os << strA.substr(u, 80) << "\n";
+ os << "\n";
+ os << ">" << nameB << "\n";
+ for(size_t u = 0; u < sz; u += 80)
+ os << strB.substr(u, 80) << "\n";
+ os << "\n";
+ } else if (format == 2) {
+ os << " " << 2 << " " << sz << "\n";
+ os << std::setw(9) << std::setiosflags(std::ios::left) << nameA << " "
+ << strA.substr(0,68) << "\n";
+
+ os << std::setw(9) << std::setiosflags(std::ios::left) << nameB << " "
+ << strB.substr(0,68) << "\n\n";
+
+ for(size_t u=68;u < sz; u+= 78) {
+ os << strA.substr(u,78) << "\n"
+ << strB.substr(u,78) << "\n\n";
+ }
+
+ } else {
+ os << "Invalid output format. Valid formats are 'aln' and 'fas'.\n";
+ }
+ os.flush();
+}
+
+#endif
diff --git a/src/config.h.cmake b/src/config.h.cmake
new file mode 100644
index 0000000..cab75b8
--- /dev/null
+++ b/src/config.h.cmake
@@ -0,0 +1,10 @@
+#pragma once
+#ifndef NGILA_CONFIG_H
+#define NGILA_CONFIG_H
+
+#cmakedefine HAVE_UNISTD_H 1
+
+#cmakedefine HAVE_GETPID 1
+
+#endif
+
diff --git a/src/lgmod.incl b/src/lgmod.incl
new file mode 100644
index 0000000..0853ef3
--- /dev/null
+++ b/src/lgmod.incl
@@ -0,0 +1,110 @@
+ 2.81186770670314445e-01, 2.36518498219483869e-01, 2.04882893380584602e-01, 2.30330197759651134e-01,
+ 1.13740933704625444e-01, 2.01908395070635932e-01, 2.67555601698039547e-01, 2.39451456458297629e-01,
+ 1.49515885443654434e-01, 2.49313056216476564e-01, 3.14771345582789686e-01, 2.54165300542776695e-01,
+ 1.51495874531288793e-01, 2.05674500120943515e-01, 2.09857094233194791e-01, 2.47380274072125650e-01,
+ 2.30839771270030497e-01, 1.09845345827668098e-01, 1.84810713975137264e-01, 2.62958171578675237e-01,
+ 5.73939347590331558e-17, -2.67181227852833103e-01, -4.05532800937279669e-01, -4.32766242741772722e-01,
+-4.87679314891660609e-01, -5.46173758272108900e-01, -6.08255144175605000e-01, -7.53916572568585464e-01,
+-9.16826265967801257e-01, -9.85658030071095981e-01, -1.00484550319211419e+00, -1.06338020029818536e+00,
+-1.16839769895865286e+00, -1.25174685113231576e+00, -1.33707818531830824e+00, -1.43352201488587894e+00,
+-1.45873407096179819e+00, -1.65342896766215786e+00, -1.79911825916068180e+00, -2.11616459749516928e+00,
+-2.81186630077034649e-01, 7.02139176999369563e-02, -9.42586187996573566e-02, 4.95215293965131729e-03,
+-1.79977829716740506e-01, -3.67515721518490504e-02, 1.18864621160894319e-02, 3.64117171017498964e-01,
+ 1.31374253388901963e-01, -3.03604759547569214e-01, -1.26268562650448057e-01, 1.27134391599006386e-01,
+ 1.09219283993825786e-01, 4.86410742317634481e-01, 5.03091279504468858e-01, 5.07475428804172249e-02,
+ 5.05637409270126362e-02, 5.70769592157559999e-02, -2.61027745320369542e-01, -1.44697133912360443e-01,
+-2.36518379960323322e-01, 1.74342611139448600e-01, 1.11788091343182613e-01, -9.41892684397988428e-02,
+ 2.83128520084429958e-01, -9.49802010565075211e-02, 5.40586312303613648e-01, -1.18871177000215461e-01,
+-8.78349480776494385e-02, 4.37380203032346898e-02, 2.28690576251423111e-03, -2.20608508020384314e-01,
+ 8.21466662407178555e-02, 4.71193340617187539e-01, -3.87677224225155215e-01, 1.54625385082458711e-01,
+ 1.84638084404970926e-01, 2.27134389189858135e-02, -2.92322083115783560e-02, -3.42626183704825005e-03,
+-2.04882790939214798e-01, 1.63492161630914346e-01, 6.73909232455021934e-02, -6.94955506604025691e-02,
+ 1.27329873525226858e-02, 3.92710293436903876e-02, -9.18307231799749268e-02, 1.01864617610535196e-01,
+ 2.02894907223298115e-01, 4.65841337730619254e-01, -1.15169809393893222e-01, -1.16153416900409728e-01,
+-1.09402247637577443e-01, -3.18013342130180454e-01, 1.63661463393964868e-01, 3.83147629048965654e-01,
+ 5.32736826214430526e-01, 7.39900590313703616e-02, -2.00966980322749678e-01, -1.09117283022541756e-02,
+-2.30330082594638552e-01, 2.73051731008476462e-01, 1.02652832030561469e-01, -1.07961342289679532e-01,
+ 8.08820106721697130e-02, 2.38380424706745031e-02, -5.48559928569121391e-01, -1.15773850857816885e-01,
+-2.94821442086192903e-02, 3.10416593493816284e-01, -5.09103245245875677e-02, -4.11830243126780571e-01,
+ 1.35517319616029669e-01, 2.60878588625666763e-01, 3.97165218524809238e-02, -7.33192293605703094e-02,
+-4.06715524320417077e-01, -2.22859610732898429e-02, 7.13240818702267926e-03, -1.18452790646117467e-03,
+-1.13740876834201254e-01, -3.40522270248994618e-02, -2.43435196860758984e-02, -1.23888749369131575e-02,
+-7.60213448587038210e-02, 1.33595448294391295e-02, 2.69441513221340030e-02, 5.04053633757999497e-01,
+-8.11605750837193107e-01, 2.01743327083754503e-01, 1.26624795270773083e-01, -3.71781766123548130e-03,
+-1.81169984033142384e-02, -8.08077183138695421e-02, -2.49034759206184161e-02, -4.90911180111759660e-03,
+-1.47293007762951504e-02, 3.75779996363324875e-03, -2.50687498840607148e-02, -2.42848818519101894e-02,
+-2.01908294116514064e-01, 1.36413385357427752e-01, 5.86924987957951894e-02, -6.09048839682131138e-02,
+ 1.57817780240436856e-01, -2.99152314680873194e-02, 4.97144283073223947e-02, -5.02558476330183801e-02,
+-9.34108599020300190e-03, -7.94258698630168070e-02, 7.84822111862326227e-02, 3.87736658553675118e-01,
+-2.26971861878986686e-01, -1.23394834027778844e-01, 2.29353370106877286e-03, 6.32519363204767759e-01,
+-5.11768905747554737e-01, -1.03406866342249784e-01, -1.67508433012589775e-02, -8.20175869586436956e-03,
+-2.67555467920339018e-01, 2.85017466815434495e-01, 1.00188235068713050e-01, -1.19400849729773384e-01,
+ 1.93504868604052971e-01, -3.83558982849062508e-02, -4.35898382429133702e-01, -1.83123058096413716e-01,
+-2.12752301130465782e-01, -4.66705121637486953e-01, 1.29478579862907062e-01, 2.60140550297849438e-01,
+-2.54985224362659496e-02, -7.94060834707280838e-02, -1.63769216638684545e-01, -1.28347945603764446e-01,
+ 4.01649481808424358e-01, 3.55646341641650404e-02, 8.70228671846628998e-03, 1.07439862055036528e-03,
+-2.39451336732659170e-01, 2.39584237742021422e-01, 2.38618940825112910e-02, -8.86376313050787346e-02,
+-8.16018508981044999e-01, 1.53673205594599738e-01, 1.87422199437619980e-01, -3.50107679119599946e-01,
+-8.21360989779986728e-02, -5.02807406227157028e-05, 9.37889001569422454e-02, -8.37149510654303097e-03,
+-2.55704966671493061e-02, -7.07555543950137422e-02, -7.64492130743325005e-02, -2.91032490697954449e-02,
+-3.43051718254586691e-02, 8.31892901541786411e-04, -2.17598461796298936e-02, 4.10035564162072180e-03,
+-1.49515810685767775e-01, 3.59829842517823129e-02, 1.16965307383584877e-01, -1.26234737953581054e-02,
+ 9.67845827956361204e-02, 1.44510837007971155e-01, 6.27355814869766520e-02, 4.99640407811032805e-02,
+ 1.76843230550260855e-01, 4.13551573606321221e-01, 2.95623775137141720e-01, 4.21837646558613222e-01,
+-4.50496819393878478e-01, 2.38856988008023990e-01, 3.47740349751114453e-02, -4.43081081982811598e-01,
+-4.55562621866201223e-03, 1.06183477520509495e-02, 2.29042047860518114e-02, 3.71319153456787371e-04,
+-2.49312931560041806e-01, -3.11261566383156774e-01, -2.42538831175779956e-01, -9.90130193917874435e-02,
+-2.03163881702027316e-03, -1.68395782098853725e-01, -4.44678519236720027e-02, -4.98608457501498736e-02,
+ 8.81576558764247953e-02, -4.17424875159563002e-02, 3.18147365756962675e-01, -2.35956879733286662e-01,
+-1.51787081666750867e-01, -7.10479290457243229e-02, -3.18678092688356451e-02, 3.62381898585496173e-02,
+ 2.29321191286557848e-02, -2.30830956266092624e-02, 1.91239545554508228e-01, -7.10669478603731863e-01,
+-3.14771188197234819e-01, -4.34948442017680459e-01, -2.34333849066825750e-01, -1.00566253344859191e-01,
+ 3.88034581864862282e-02, -1.26398133120823541e-01, -4.03171331152079732e-02, -3.24631234597722751e-01,
+-1.28106518012837844e-01, 2.69792700925987283e-01, -2.64533985350441092e-01, 3.52137889251615821e-01,
+ 3.52289817413702910e-01, 3.60613277509554686e-02, 2.74853804406818003e-02, -5.24571255144215820e-02,
+ 6.01310341768210022e-02, -3.02386306098948399e-01, -8.42399435572435407e-02, 8.93078903925503542e-02,
+-2.54165173460221661e-01, 2.10676096678484026e-01, 9.55528152970522293e-02, -1.07194014247812427e-01,
+ 2.80265358835397771e-01, -9.65988323537882998e-02, 3.96168589250013359e-01, -1.08497155656001265e-01,
+-5.63237626652494106e-02, -9.09301714163790881e-02, -5.27892886697688343e-02, -1.55651244076054962e-01,
+ 7.73793723611961026e-02, -4.33729356251832510e-01, 4.75737025405233549e-01, -3.58093446172875862e-01,
+-1.55786644258408158e-01, -2.64511830180615415e-02, 4.10223648566811483e-02, 4.64172558238677094e-03,
+-1.51495798783408298e-01, -1.54658466299068320e-01, -9.28023859385347649e-02, -4.87790768548515830e-02,
+ 1.53682871074888039e-02, -6.18640153004653692e-02, -3.40048228774950112e-03, -7.30176028824163215e-02,
+-1.13175673207028815e-02, 6.57855605866906890e-02, -7.15552939875762628e-02, 1.23186075802851069e-01,
+ 1.06745191100497988e-01, -3.14209620082874508e-02, -2.45366492326655750e-02, 1.94548373800332297e-02,
+-1.07479497544866853e-01, 9.25357775500510460e-01, 1.38435988470992505e-01, 8.72645468883845865e-02,
+-2.05674397283770571e-01, -3.47882764068484396e-01, 1.81761524912582545e-01, 4.43553193046316852e-02,
+ 4.55608697460548329e-02, 4.18243628872892292e-01, -1.95551856991488249e-02, -9.54377701379390264e-02,
+-1.39129663819227561e-01, -1.95875241924726912e-01, -5.09592905960478459e-01, -2.03137357832441967e-01,
+-5.04503312545028781e-01, 3.81002591725322082e-02, -4.69428441680909414e-03, -1.80483023700174702e-02,
+-7.24441374690579232e-03, 8.32361890322731897e-03, -7.58000998451296894e-03, -7.54462355819368559e-04,
+-2.09856989304726321e-01, 1.69761411548488478e-01, -3.15231110987661378e-01, 8.91175006278178694e-01,
+ 8.39474034380392892e-02, 6.12102334722513519e-02, 1.08333972461327477e-03, -1.24020099203862477e-01,
+-4.71639536951998783e-02, 3.51945829518460762e-02, 3.59877092487167344e-02, -3.32415241116488652e-02,
+-1.02924402059669795e-02, -2.82875220101856950e-02, -1.50334201307824669e-02, 2.15526400278748254e-04,
+ 6.36546052612529945e-03, 9.69173937711236520e-03, -1.62170358743539862e-02, -2.95352181095841866e-03,
+-2.47380150382081249e-01, 1.21459608399252461e-01, -2.57137128478748464e-02, 4.58813916137140732e-03,
+-1.02099238910233522e-01, 1.09699591463227608e-03, 4.48135996086101995e-03, 3.47709883859115221e-01,
+ 1.91241457150310490e-01, -1.07840601453720875e-03, -2.48270760076264779e-01, 7.17038586356512075e-02,
+ 8.70247260410817319e-02, -2.58552329878836938e-02, -1.18047846110177546e-01, 7.25130165942866611e-03,
+ 5.41569149096923219e-02, -1.56753837739141477e-01, 7.90977794760870379e-01, 1.15265445725544799e-01,
+-2.30839655850231390e-01, 4.16813430839604454e-02, -6.58207161847328415e-02, -2.67337745228739543e-02,
+-3.00659792705283771e-02, -5.67035933465922792e-02, 9.77688869888826417e-03, 3.58549985085491685e-01,
+ 2.75051780197904661e-01, -6.40816365884891576e-02, -2.00018328698603781e-01, 2.52840455233739603e-02,
+ 6.52568835988332485e-02, -2.78710843542326892e-01, -5.38702355822758605e-01, -2.57576512718465822e-01,
+-2.15375790270600553e-01, -1.62790120009007793e-03, -4.41363996644431045e-01, -8.41042388970196508e-02,
+-1.09845290905036308e-01, -2.31953391524271413e-01, 7.00047946640009910e-01, 3.20931180945844496e-01,
+-1.70058541102086014e-01, -5.55571187437136516e-01, -5.38101766240506874e-02, 6.74942482949529750e-03,
+ 1.31283582215044672e-02, 1.01428560755042751e-02, 9.13168551794209943e-03, 2.03339139592957630e-03,
+-8.06495675145409152e-03, -7.05015252981814472e-03, 6.66338213154941093e-03, -3.75628077872242061e-03,
+ 9.22440873091096014e-05, -2.47358435492690170e-03, -6.53779112287955467e-04, -2.18573828768679887e-03,
+-1.84810621569849537e-01, -2.52693577559623628e-01, 3.45189755101757312e-01, 8.62435990434276406e-02,
+ 8.19894764354378996e-02, 6.02564606926505042e-01, 2.04664398108178602e-02, 8.90973282552566942e-02,
+ 1.24722762340539972e-01, -6.09510242181348708e-02, 3.82610409617307523e-01, -6.87604247279448300e-03,
+ 4.67680500158018653e-01, -7.70560724311262074e-02, -3.20336558848138137e-03, 8.18597282502086243e-02,
+-6.80203127957003384e-03, -9.43472150911189161e-03, -3.78781691994870516e-03, -1.03977614870460621e-03,
+-2.62958040099687873e-01, -2.55605372713718859e-01, -2.38976150750236899e-01, -8.58208808017060509e-02,
+-3.11584620773551063e-02, -1.71311742071191597e-01, -4.03595556073759187e-02, 8.60147616421999550e-02,
+ 1.17880542571498881e-01, -1.48456189850605325e-01, 3.80584544724422991e-01, -2.97906786753040387e-01,
+-2.12397023579544542e-01, -5.17215199470540642e-03, 2.46282337296211684e-02, 3.84507666838677375e-02,
+ 7.73619319858645516e-03, -6.61815380122089819e-02, -5.73496046298600237e-02, 6.61223828104091194e-01
diff --git a/src/matparser.cpp b/src/matparser.cpp
new file mode 100644
index 0000000..cf18124
--- /dev/null
+++ b/src/matparser.cpp
@@ -0,0 +1,117 @@
+/****************************************************************************
+ * Copyright (C) 2005-2010 Reed A. Cartwright, PhD <reed at scit.us> *
+ * *
+ * This program is free software: you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation, either version 3 of the License, or *
+ * (at your option) any later version. *
+ * *
+ * This program is distributed in the hope that it will be useful, *
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of *
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
+ * GNU General Public License for more details. *
+ * *
+ * You should have received a copy of the GNU General Public License *
+ * along with this program. If not, see <http://www.gnu.org/licenses/>. *
+ ****************************************************************************/
+
+#include "ngila.h"
+
+#include "matparser.h"
+
+using namespace std;
+
+bool parse_matrix(const char *cs, sub_matrix &rsm, bool bi)
+{
+ mat_work w;
+ mat_grammar my_grammar(w);
+ if(strcmp(cs, "-")==0)
+ {
+ string ss;
+ getline(cin, ss, '\004');
+ if(ss.empty())
+ return CERRORR("unable to open stdin.");
+ parse_info<string::const_iterator> info = parse(ss.begin(), ss.end(), my_grammar, blank_p);
+ if (!info.full)
+ return CERRORR("unable to parse stdin.");
+ }
+ else
+ {
+ file_iterator<char> file_first(cs);
+ if(!file_first)
+ return CERRORR("unable to open \'" << cs << "\'.");
+ file_iterator<char> file_last = file_first.make_end();
+ parse_info< file_iterator<char> > info = parse(file_first, file_last, my_grammar, blank_p);
+ if (!info.full)
+ return CERRORR("unable to parse \'" << cs << "\'.");
+ }
+ if(!w.process(rsm, bi))
+ return CERRORR("matrix \'" << cs << "\' must be triangular or square.");
+ return true;
+}
+
+bool mat_work::process(matrix &m, bool bi) const
+{
+ sub_matrix_clear(m);
+
+ size_t sz = labels.size();
+
+ if(data.size() != sz || !(data.front().size() == sz
+ || data.back().size() == sz) )
+ return false;
+ labels_type llabels(labels);
+ if(bi)
+ {
+ for(labels_type::iterator it = llabels.begin();
+ it != llabels.end(); ++it)
+ *it = toupper(*it);
+ }
+ if(data.front().size() == 1)
+ {
+ // Assume lower triangular matrix
+ for(size_t i = 0; i < sz; ++i)
+ {
+ if(data[i].size() != i+1)
+ return false;
+ size_t ich = static_cast<size_t>(llabels[i]);
+ for(size_t j = 0; j <= i; ++j)
+ {
+ size_t jch = static_cast<size_t>(llabels[j]);
+ m[ich][jch] = m[jch][ich] = data[i][j];
+ }
+ }
+ }
+ else if(data.back().size() == 1)
+ {
+ // Assume upper triangular matrix
+ for(size_t i = 0; i < sz; ++i)
+ {
+ if(data[i].size() != sz-i)
+ return false;
+ size_t ich = static_cast<size_t>(llabels[i]);
+ for(size_t j = 0; j < sz-i; ++j)
+ {
+ size_t jch = static_cast<size_t>(llabels[j+i]);
+ m[ich][jch] = m[jch][ich] = data[i][j];
+ }
+ }
+ }
+ else
+ {
+ // Assume square matrix
+ for(size_t i = 0; i < sz; ++i)
+ {
+ if(data[i].size() != sz)
+ return false;
+ size_t ich = static_cast<size_t>(llabels[i]);
+ for(size_t j = 0; j < sz; ++j)
+ {
+ size_t jch = static_cast<size_t>(llabels[j]);
+ m[ich][jch] = data[i][j];
+ }
+ }
+ }
+ return true;
+}
+
+
diff --git a/src/matparser.h b/src/matparser.h
new file mode 100644
index 0000000..89d775a
--- /dev/null
+++ b/src/matparser.h
@@ -0,0 +1,92 @@
+/****************************************************************************
+ * Copyright (C) 2005-2010 Reed A. Cartwright, PhD <reed at scit.us> *
+ * *
+ * This program is free software: you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation, either version 3 of the License, or *
+ * (at your option) any later version. *
+ * *
+ * This program is distributed in the hope that it will be useful, *
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of *
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
+ * GNU General Public License for more details. *
+ * *
+ * You should have received a copy of the GNU General Public License *
+ * along with this program. If not, see <http://www.gnu.org/licenses/>. *
+ ****************************************************************************/
+
+#ifndef MATPARSER_H
+#define MATPARSER_H
+
+#include <limits>
+
+#include <boost/spirit/include/classic.hpp>
+#include <boost/spirit/include/classic_core.hpp>
+#include <boost/spirit/include/classic_parser.hpp>
+#include <boost/spirit/include/classic_assign_actor.hpp>
+#include <boost/spirit/include/classic_clear_actor.hpp>
+using namespace boost::spirit::classic;
+
+typedef double sub_matrix[128][128];
+const size_t sub_matrix_size = 128;
+inline void sub_matrix_clear(sub_matrix &m) {
+ std::fill(&m[0][0], (&m[0][0])+sub_matrix_size*sub_matrix_size,
+ std::numeric_limits<double>::max()/16.0);
+ for(size_t i=0;i<sub_matrix_size;++i)
+ m[i][i] = 0.0;
+}
+
+bool parse_matrix(const char *cs, sub_matrix &rsm, bool bi = false);
+
+struct mat_work
+{
+ typedef std::vector<char> labels_type;
+ typedef std::vector<double> row_type;
+ typedef std::vector<row_type> data_type;
+ typedef sub_matrix matrix;
+
+ labels_type labels;
+ data_type data;
+ row_type rt;
+
+ inline void clear()
+ {
+ labels.clear();
+ data.clear();
+ rt.clear();
+ }
+
+ bool process(matrix &m, bool bi = false) const;
+};
+
+struct mat_grammar : public grammar<mat_grammar>
+{
+ mat_grammar(mat_work &w) : work(w) { work.clear(); }
+
+ mat_work &work;
+
+ template <typename ScannerT> struct definition
+ {
+ definition(mat_grammar const& self)
+ {
+ matrix = *eol_p >> header >> body >> end_p;
+ header = +(graph_p[push_back_a(self.work.labels)]);
+ body = +line >> !eol_p;
+ line = eol_p >> !(graph_p - real_p) >>
+ (+(real_p[push_back_a(self.work.rt)]))
+ [push_back_a(self.work.data, self.work.rt)]
+ [clear_a(self.work.rt)]
+ ;
+ }
+
+ rule<ScannerT> matrix;
+ rule<ScannerT> header;
+ rule<ScannerT> body;
+ rule<ScannerT> line;
+ rule<ScannerT> const& start() const { return matrix; }
+ };
+};
+
+#endif
+
+
diff --git a/src/models.cpp b/src/models.cpp
new file mode 100644
index 0000000..d14906a
--- /dev/null
+++ b/src/models.cpp
@@ -0,0 +1,442 @@
+/****************************************************************************
+ * Copyright (C) 2005-2010 Reed A. Cartwright, PhD <reed at scit.us> *
+ * *
+ * This program is free software: you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation, either version 3 of the License, or *
+ * (at your option) any later version. *
+ * *
+ * This program is distributed in the hope that it will be useful, *
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of *
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
+ * GNU General Public License for more details. *
+ * *
+ * You should have received a copy of the GNU General Public License *
+ * along with this program. If not, see <http://www.gnu.org/licenses/>. *
+ ****************************************************************************/
+
+#include "ngila.h"
+#include "models.h"
+
+#include <limits>
+#include <iomanip>
+
+#ifdef USE_GSL_ZETA
+#include <gsl/gsl_sf_zeta.h>
+inline double log_zeta(double z) { return log(gsl_sf_zeta(z)); }
+#else
+double log_zeta(double z) {
+ // spline interpolation of
+ // y = log((z-1)*zeta(z)) for z = seq(1,10.1,0.1)
+ static double y[92] = {
+ 0.00000000000000000e+00, 5.68007048928620834e-02, 1.11824418961510755e-01, 1.65162481245235387e-01,
+ 2.16899224977457172e-01, 2.67112722170841699e-01, 3.15875427900892403e-01, 3.63254740373157525e-01,
+ 4.09313489900198446e-01, 4.54110367565164985e-01, 4.97700302470745981e-01, 5.40134794961248610e-01,
+ 5.81462211982759714e-01, 6.21728049750888800e-01, 6.60975168080968012e-01, 6.99244000065426419e-01,
+ 7.36572740229037959e-01, 7.72997513832677563e-01, 8.08552529612340476e-01, 8.43270217918634102e-01,
+ 8.77181355951436870e-01, 9.10315181555988784e-01, 9.42699496853069863e-01, 9.74360762811215531e-01,
+ 1.00532418572836546e+00, 1.03561379646993168e+00, 1.06525252320689856e+00, 1.09426225830844959e+00,
+ 1.12266391996660997e+00, 1.15047750906370894e+00, 1.17772216173544519e+00, 1.20441619803187527e+00,
+ 1.23057716703450559e+00, 1.25622188874908414e+00, 1.28136649305977168e+00, 1.30602645600060807e+00,
+ 1.33021663357389319e+00, 1.35395129332193020e+00, 1.37724414383803273e+00, 1.40010836238450520e+00,
+ 1.42255662076911871e+00, 1.44460110961722821e+00, 1.46625356116381989e+00, 1.48752527067832974e+00,
+ 1.50842711662481577e+00, 1.52896957965086511e+00, 1.54916276049037749e+00, 1.56901639685795402e+00,
+ 1.58853987940593444e+00, 1.60774226680911636e+00, 1.62663230003675863e+00, 1.64521841586655548e+00,
+ 1.66350875969082224e+00, 1.68151119766111612e+00, 1.69923332821384077e+00, 1.71668249301607734e+00,
+ 1.73386578736785757e+00, 1.75079007009433818e+00, 1.76746197295883389e+00, 1.78388790962537036e+00,
+ 1.80007408419733039e+00, 1.81602649935684246e+00, 1.83175096412782179e+00, 1.84725310128395459e+00,
+ 1.86253835442144466e+00, 1.87761199471499252e+00, 1.89247912737422341e+00, 1.90714469781663887e+00,
+ 1.92161349757210820e+00, 1.93589016993294716e+00, 1.94997921536272623e+00, 1.96388499667612382e+00,
+ 1.97761174400136919e+00, 1.99116355953612101e+00, 2.00454442210694284e+00, 2.01775819154196467e+00,
+ 2.03080861286571279e+00, 2.04369932032459678e+00, 2.05643384125103657e+00, 2.06901559977375671e+00,
+ 2.08144792038136428e+00, 2.09373403134591785e+00, 2.10587706801283314e+00, 2.11788007596312688e+00,
+ 2.12974601405368036e+00, 2.14147775734089807e+00, 2.15307809989286270e+00, 2.16454975749481671e+00,
+ 2.17589537025256385e+00, 2.18711750509813818e+00, 2.19821865820188833e+00, 2.20920125729490113e+00
+ };
+ static double b[92] = {
+ 5.77198067017487992e-01, 5.58968262571901242e-01, 5.41661451540230487e-01, 5.25239221838376791e-01,
+ 5.09625841584655737e-01, 4.94764639591190325e-01, 4.80601687753637097e-01, 4.67089155463733474e-01,
+ 4.54183550370610756e-01, 4.41845458814050907e-01, 4.30038991489612388e-01, 4.18731397110005898e-01,
+ 4.07892705430773550e-01, 3.97495424856108670e-01, 3.87514278091038833e-01, 3.77925972215862271e-01,
+ 3.68708997487613543e-01, 3.59843450851215829e-01, 3.51310880606601483e-01, 3.43094149301072793e-01,
+ 3.35177312361997326e-01, 3.27545510371581250e-01, 3.20184873200665854e-01, 3.13082434482556182e-01,
+ 3.06226055127975738e-01, 2.99604354767027625e-01, 2.93206650159908988e-01, 2.87022899748871851e-01,
+ 2.81043653635944191e-01, 2.75260008365130482e-01, 2.69663565968592922e-01, 2.64246396805493722e-01,
+ 2.59001005781242910e-01, 2.53920301585792152e-01, 2.48997568633569616e-01, 2.44226441425652752e-01,
+ 2.39600881087469908e-01, 2.35115153864130427e-01, 2.30763811380187006e-01, 2.26541672492370061e-01,
+ 2.22443806582916848e-01, 2.18465518157651667e-01, 2.14602332627510634e-01, 2.10849983165356192e-01,
+ 2.07204398540939916e-01, 2.03661691846944048e-01, 2.00218150038134451e-01, 1.96870224213184369e-01,
+ 1.93614520575835614e-01, 1.90447792018342332e-01, 1.87366930275524757e-01, 1.84368958602731337e-01,
+ 1.81451024935457061e-01, 1.78610395492258517e-01, 1.75844448786063728e-01, 1.73150670012326841e-01,
+ 1.70526645785132130e-01, 1.67970059194969018e-01, 1.65478685164280659e-01, 1.63050386078872789e-01,
+ 1.60683107675126668e-01, 1.58374875164782741e-01, 1.56123789580483313e-01, 1.53928024326647012e-01,
+ 1.51785821921613967e-01, 1.49695490918038104e-01, 1.47655402989595230e-01, 1.45663990172970714e-01,
+ 1.43719742255064759e-01, 1.41821204296018161e-01, 1.39966974279406492e-01, 1.38155700881651272e-01,
+ 1.36386081353280036e-01, 1.34656859505147358e-01, 1.32966823793335631e-01, 1.31314805496822951e-01,
+ 1.29699676982467121e-01, 1.28120350052274373e-01, 1.26575774368151456e-01, 1.25064935949913525e-01,
+ 1.23586855742028531e-01, 1.22140588246809095e-01, 1.20725220214797091e-01, 1.19339869410269847e-01,
+ 1.17983683369542569e-01, 1.16655838444697929e-01, 1.15355538027138341e-01, 1.14082014064304196e-01,
+ 1.12834516506675933e-01, 1.11612348008638756e-01, 1.10414729938505748e-01, 1.09241298140229107e-01
+ };
+ static double c[92] = {
+ -9.34324982101620849e-02, -8.88655462457056400e-02, -8.42025640710015594e-02, -8.00197329475356850e-02,
+ -7.61140695896742742e-02, -7.24979503449798279e-02, -6.91315680305528008e-02, -6.59937548684829400e-02,
+ -6.30622960627443691e-02, -6.03186195028546090e-02, -5.77460537415301622e-02, -5.53298900545345917e-02,
+ -5.30570267377887678e-02, -5.09157790088607815e-02, -4.88956886418370762e-02, -4.69873701099287749e-02,
+ -4.51823771725581441e-02, -4.34730891914191925e-02, -4.18526132547244814e-02, -4.03146998005620635e-02,
+ -3.88536695901928741e-02, -3.74643503139678175e-02, -3.61420213951862412e-02, -3.48823657859102590e-02,
+ -3.36814277598939601e-02, -3.25355758495872241e-02, -3.14414702215991201e-02, -3.03960338887724092e-02,
+ -2.93964272405043761e-02, -2.84400254676326451e-02, -2.75243984977431091e-02, -2.66472931332491295e-02,
+ -2.58066171092583355e-02, -2.50004248452490169e-02, -2.42269046769766766e-02, -2.34843674021920808e-02,
+ -2.27712359796366053e-02, -2.20860362537580353e-02, -2.14273885856757704e-02, -2.07940002924938278e-02,
+ -2.01846588020384680e-02, -1.95982254506131687e-02, -1.90336298507971825e-02, -1.84898647707473392e-02,
+ -1.79659814734153059e-02, -1.74610854665435272e-02, -1.69743326215522836e-02, -1.65049256279485614e-02,
+ -1.60521107455388889e-02, -1.56151748293941583e-02, -1.51934425987815744e-02, -1.47862741291524755e-02,
+ -1.43930625435904746e-02, -1.40132318883947268e-02, -1.36462351735532200e-02, -1.32915525638158481e-02,
+ -1.29486897081312058e-02, -1.26171761934997899e-02, -1.22965641133836890e-02, -1.19864267406951872e-02,
+ -1.16863572967661634e-02, -1.13959678066728187e-02, -1.11148880363215272e-02, -1.08427645020414097e-02,
+ -1.05792595482891699e-02, -1.03240504874694411e-02, -1.00768287969594007e-02, -9.83729936928581100e-03,
+ -9.60517980977360075e-03, -9.38019978069253024e-03, -9.16210038542411738e-03, -8.95063359212771269e-03,
+ -8.74556169158521510e-03, -8.54665678974113850e-03, -8.35370032837635598e-03, -8.16648263675044672e-03,
+ -7.98480250680756731e-03, -7.80846679512038448e-03, -7.63729004610835815e-03, -7.47109413627102155e-03,
+ -7.30970794257910974e-03, -7.15296700961526831e-03, -7.00071331050457386e-03, -6.85279473476791898e-03,
+ -6.70906567250490436e-03, -6.56938357594156162e-03, -6.43362059965430465e-03, -6.30161902868694133e-03,
+ -6.17335654759581381e-03, -6.04832843277599772e-03, -5.92785226855403573e-03, -5.80646571421250274e-03
+ };
+ static double d[92] = {
+ 1.52231732148548141e-02, 1.55432739156802779e-02, 1.39427704115529256e-02, 1.30188778595380535e-02,
+ 1.20537308156481206e-02, 1.12212743814234135e-02, 1.04593772068995155e-02, 9.77152935246192872e-03,
+ 9.14558853296590195e-03, 8.57521920441479263e-03, 8.05387895665191783e-03, 7.57621105581938320e-03,
+ 7.13749242975999110e-03, 6.73363455674563707e-03, 6.36106177302769339e-03, 6.01664312456877577e-03,
+ 5.69762660379651183e-03, 5.40158645564904442e-03, 5.12637818054137179e-03, 4.87010070123063132e-03,
+ 4.63106425408353530e-03, 4.40776306260523814e-03, 4.19885203091994313e-03, 4.00312675338765356e-03,
+ 3.81950636768913941e-03, 3.64701875996034349e-03, 3.48478777608903336e-03, 3.33202216089343488e-03,
+ 3.18800590957243933e-03, 3.05208989963179746e-03, 2.92368454831327604e-03, 2.80225341330263191e-03,
+ 2.68730754669771433e-03, 2.57840056090781001e-03, 2.47512424928199477e-03, 2.37710474185159329e-03,
+ 2.28399908626188489e-03, 2.19549222694087999e-03, 2.11129431060647988e-03, 2.03113830151787017e-03,
+ 1.95477783808431504e-03, 1.88198533271996342e-03, 1.81255026683282061e-03, 1.74627765777343497e-03,
+ 1.68298668957260738e-03, 1.62250948330413397e-03, 1.56468997867907677e-03, 1.50938294136556956e-03,
+ 1.45645305381577120e-03, 1.40577410204195132e-03, 1.35722823209698923e-03, 1.31070528520667415e-03,
+ 1.26610218398581624e-03, 1.22332238280502704e-03, 1.18227536579124689e-03, 1.14287618561547104e-03,
+ 1.10504504877138735e-03, 1.06870693372033065e-03, 1.03379124229500984e-03, 1.00023147976341462e-03,
+ 9.67964966977814638e-04, 9.36932567837640211e-04, 9.07078447600385431e-04, 8.78349845840801211e-04,
+ 8.50696869399103234e-04, 8.24072301700124322e-04, 7.98431425578638139e-04, 7.73731865040699603e-04,
+ 7.49933430270234851e-04, 7.26997984228045478e-04, 7.04889310988004222e-04, 6.83573001808341274e-04,
+ 6.63016339480245885e-04, 6.43188204549280233e-04, 6.24058972086363479e-04, 6.05600433142925690e-04,
+ 5.87785705623952454e-04, 5.70589163373411514e-04, 5.53986366124460195e-04, 5.37953978973039872e-04,
+ 5.22469776546138497e-04, 5.07512330368975668e-04, 4.93061919122183244e-04, 4.79096874210053288e-04,
+ 4.65606988544477478e-04, 4.52543254290855321e-04, 4.40005236557871476e-04, 4.27541603637093251e-04,
+ 4.16760382732723241e-04, 4.01587214073209479e-04, 4.04621847805109239e-04, 4.04621847805109239e-04
+ };
+ // truncate at z=11
+ if(z > 11.0)
+ z = 11.0;
+ // if z is out of range, just return 0.0
+ if(z <= 1.0)
+ return 0.0;
+ // find position in table
+ size_t p = static_cast<size_t>((z-1.0)*10.0);
+ if(p > 91)
+ p = 91;
+ double dx = z-(1.0+p*0.1);
+ // spline approximation
+ return y[p] + dx*(b[p] + dx*(c[p] + dx*d[p]))
+ - log(z-1.0);
+}
+#endif
+
+using namespace std;
+
+/****************************************************************************
+ * class cost_model *
+ ****************************************************************************/
+
+
+bool cost_model::create(const ngila_app::args &rargs)
+{
+ dA = rargs.cost_intersection;
+ dB = rargs.cost_linear;
+ dC = rargs.cost_logarithmic;
+ dF = rargs.cost_intersection_free;
+ dG = rargs.cost_linear_free;
+ dH = rargs.cost_logarithmic_free;
+
+ if(rargs.cost_matrix.empty())
+ {
+ for(size_t i=0;i<sub_matrix_size;++i)
+ for(size_t j=0;j<sub_matrix_size;++j)
+ mCost[i][j] = (i == j) ? rargs.cost_match : rargs.cost_mismatch;
+ }
+ else
+ {
+ if(!parse_matrix(rargs.cost_matrix.c_str(), mCost, rargs.case_insensitivity))
+ return CERRORR("parsing of \'" << rargs.cost_matrix.c_str() << "\' failed.");
+ }
+ return true;
+}
+
+/****************************************************************************
+ * class k2p_model *
+ ****************************************************************************/
+
+enum {nA, nC, nG, nT, nN, nSize};
+const int nuc_table[] = {
+ /*64*/ -1, nA, nN, nC, nN, -1, -1, nG, nN, -1, -1, nN, -1, nN, nN, -1,
+ /*80*/ -1, -1, nN, nN, nT, nT, nN, nN, -1, nN, -1, -1, -1, -1, -1, -1,
+ /*96*/ -1, nA, nN, nC, nN, -1, -1, nG, nN, -1, -1, nN, -1, nN, nN, -1,
+ /*112*/ -1, -1, nN, nN, nT, nT, nN, nN, -1, nN, -1, -1, -1, -1, -1, -1,
+ /*128*/
+};
+const int pupy[] = {1, 0, 1, 0, -1};
+
+bool k2p_model::create(const ngila_app::args &rargs)
+{
+ if(rargs.branch_length <= 0.0)
+ return CERRORR("branch length must be positive.");
+ if(rargs.ratio <= 0.0)
+ return CERRORR("Ts/Tv ratio must be positive.");
+ if(rargs.avgaln <= 0.0)
+ return CERRORR("avgaln must be positive.");
+ if(rargs.indel_rate <= 0.0)
+ return CERRORR("indel rate must be positive.");
+
+ double p_ts = 0.25-0.5*exp(-rargs.branch_length*(2.0*rargs.ratio+1.0)/(rargs.ratio+1.0))
+ + 0.25*exp(-2.0*rargs.branch_length/(rargs.ratio+1.0));
+ double p_tv = 0.5-0.5*exp(-2.0*rargs.branch_length/(rargs.ratio+1.0));
+ double p_match = 1.0-p_ts-p_tv;
+ double l_h = -2.0*rargs.indel_rate*rargs.branch_length
+ +log(rargs.avgaln)-log(rargs.avgaln+1.0);
+
+ double c_ts = -(log(p_ts)+l_h+log(0.25));
+ double c_tv = -(log(p_tv)+l_h+log(0.125));
+ double c_m = -(log(p_match)+l_h+log(0.25));
+
+ dNucScale = rargs.no_scale ? 0.0 : 0.5*c_m;
+
+ double sub_costs[nSize][nSize];
+ for(size_t i = 0;i<nN;++i)
+ {
+ sub_costs[i][i] = c_m;
+ for(size_t j = 0; j < i; ++j)
+ {
+ sub_costs[i][j] = sub_costs[j][i] = (pupy[i] == pupy[j]) ? c_ts : c_tv;
+ }
+ sub_costs[nN][i] = sub_costs[i][nN] = -(l_h+log(0.0625)); // All N's are weighted by 0.25
+ }
+ sub_costs[nN][nN] = -(l_h+log(0.0625)); // All N's are weighted by 0.25
+
+ sub_matrix_clear(mCost);
+ for(size_t i = 0;i<64;++i)
+ {
+ int ni = nuc_table[i];
+ for(size_t j=0;j<64;++j)
+ {
+ int nj = nuc_table[j];
+ if(ni == -1 || nj == -1)
+ continue;
+ mCost[i+64][j+64] = sub_costs[ni][nj] - 2.0*dNucScale;
+ }
+ }
+ dEnd = rargs.no_scale ? log(rargs.avgaln+1.0) : 0.0;
+
+ return true;
+}
+
+template<class T>
+struct isambnuc
+{
+ bool operator()(T n)
+ {
+ return (nuc_table[n-64] == nN);
+ }
+};
+
+double k2p_model::offset(const std::string &seqA, const std::string &seqB) const
+{
+ isambnuc<std::string::value_type> cpred;
+ size_t ncount = count_if(seqA.begin(), seqA.end(), cpred);
+ ncount += count_if(seqB.begin(), seqB.end(), cpred);
+
+ return dEnd-static_cast<double>(ncount)*log(4.0);
+}
+
+/****************************************************************************
+ * class zeta_model *
+ ****************************************************************************/
+
+bool zeta_model::create(const ngila_app::args &rargs) {
+ if(!k2p_model::create(rargs))
+ return false;
+ if(rargs.indel_slope <= 1.0)
+ return CERRORR("indel slope must be greater than 1.0");
+
+ dA = -(log(0.5)+log(1.0-exp(-2.0*rargs.indel_rate*rargs.branch_length))
+ + log(rargs.avgaln)-log(rargs.avgaln+1.0)
+ - log_zeta(rargs.indel_slope));
+ dB = -log(0.25) - dNucScale; // All N's are weighted by 0.25
+ dC = rargs.indel_slope;
+
+ dF = dH = 0.0;
+ //dG = -log(0.25) - dNucScale;
+ dG = rargs.no_scale ? -log(0.25) : 0.0;
+
+ return true;
+}
+
+/****************************************************************************
+ * class geo_model *
+ ****************************************************************************/
+
+bool geo_model::create(const ngila_app::args &rargs) {
+ if(!k2p_model::create(rargs))
+ return false;
+ if(rargs.indel_mean <= 1.0)
+ return CERRORR("indel mean must be greater than 1.0");
+
+ dA = -(log(0.5)+log(1.0-exp(-2.0*rargs.indel_rate*rargs.branch_length))
+ + log(rargs.avgaln)-log(rargs.avgaln+1.0)
+ - log(rargs.indel_mean-1));
+ dB = -(log(0.25)+log(rargs.indel_mean-1.0)-log(rargs.indel_mean)) - dNucScale; // All N's are weighted by 0.25
+ dC = 0.0;
+
+ dF = dH = 0.0;
+ //dG = -log(0.25) - dNucScale;
+ dG = rargs.no_scale ? -log(0.25) : 0.0;
+
+ return true;
+}
+
+/****************************************************************************
+ * class aa_model *
+ ****************************************************************************/
+
+bool aa_model::create(const ngila_app::args &rargs) {
+ enum {D=0,E=20,U=40};
+ static const char AA[] = "ARNDCQEGHILKMFPSTWYV";
+ static const char aa[] = "arndcqeghilkmfpstwyv";
+
+ static const double data[440] = {
+# include "lgmod.incl"
+ };
+ if(rargs.branch_length <= 0.0)
+ return CERRORR("branch length must be positive.");
+ if(rargs.avgaln <= 0.0)
+ return CERRORR("avgaln must be positive.");
+ if(rargs.indel_rate <= 0.0)
+ return CERRORR("indel rate must be positive.");
+
+ // skeleton parameters
+ double l_h = -2.0*rargs.indel_rate*rargs.branch_length
+ +log(rargs.avgaln)-log(rargs.avgaln+1.0);
+
+
+ double el[20], x[20][20], m[20];
+ sub_matrix_clear(mCost);
+
+ // cost of stationary amino acids
+ std::fill(&vAACost[0], &vAACost[128],
+ std::numeric_limits<double>::max()/16.0);
+ for(int i=0;i<20;++i)
+ vAACost[AA[i]] = vAACost[aa[i]] = -2.0*log(data[D+i]);
+
+ // exp the eigenvalues
+ for(int i=0;i<20;++i)
+ el[i] = exp(rargs.branch_length*data[E+i]);
+ // D*U*el
+ for(int i=0;i<20;++i) {
+ for(int j=0;j<20;++j) {
+ x[i][j] = data[U+j+20*i]*el[j]/data[D+i];
+ }
+ }
+ // x*(U^T*D^-1)
+ // matches
+ for(int i=0;i<20;++i) {
+ double scost = 0.0;
+ for(int k=0;k<20;++k)
+ scost += x[i][k]*data[U+k+20*i]*data[D+i];
+ scost = -(log(scost))-l_h-vAACost[AA[i]];
+ m[i] = 0.0;
+ if(!rargs.no_scale) {
+ m[i] = scost+2*vAACost[AA[i]];
+ scost = -2*vAACost[AA[i]];
+ }
+ mCost[AA[i]][AA[i]] = mCost[AA[i]][aa[i]]
+ = mCost[aa[i]][AA[i]] = mCost[aa[i]][aa[i]]
+ = scost;
+ }
+ // mismatches
+ for(int i=0;i<20;++i) {
+ for(int j=i+1;j<20;++j) {
+ double scost = 0.0;
+ for(int k=0;k<20;++k)
+ scost += x[i][k]*data[U+k+20*j]*data[D+j];
+ scost = -(log(scost))-l_h-vAACost[AA[j]];
+ scost -= (m[i]+m[j])/2;
+
+ mCost[AA[i]][AA[j]] = mCost[AA[j]][AA[i]]
+ = mCost[AA[i]][aa[j]] = mCost[AA[j]][aa[i]]
+ = mCost[aa[i]][AA[j]] = mCost[aa[j]][AA[i]]
+ = mCost[aa[i]][aa[j]] = mCost[aa[j]][aa[i]]
+ = scost;
+ }
+ }
+ dEnd = rargs.no_scale ? log(rargs.avgaln+1.0) : 0.0;
+
+ return true;
+}
+
+double aa_model::offset(const string &seqA, const string &seqB) const {
+ double off = 0.0;
+ for(string::size_type k=0;k<seqA.size();++k) {
+ off += vAACost[seqA[k]];
+ }
+ for(string::size_type k=0;k<seqB.size();++k) {
+ off += vAACost[seqB[k]];
+ }
+ return dEnd+off;
+}
+
+/****************************************************************************
+ * class aazeta_model *
+ ****************************************************************************/
+
+bool aazeta_model::create(const ngila_app::args &rargs) {
+ if(!aa_model::create(rargs))
+ return false;
+ if(rargs.indel_slope <= 1.0)
+ return CERRORR("indel slope must be greater than 1.0");
+
+ dA = -(log(0.5)+log(1.0-exp(-2.0*rargs.indel_rate*rargs.branch_length))
+ + log(rargs.avgaln)-log(rargs.avgaln+1.0)
+ - log_zeta(rargs.indel_slope));
+ dB = 0.0;
+ dC = rargs.indel_slope;
+ dF = dH = 0.0;
+ dG = 0.0;
+
+ return true;
+}
+
+/****************************************************************************
+ * class aageo_model *
+ ****************************************************************************/
+
+bool aageo_model::create(const ngila_app::args &rargs) {
+ if(!aa_model::create(rargs))
+ return false;
+ if(rargs.indel_mean <= 1.0)
+ return CERRORR("indel mean must be greater than 1.0");
+
+ dA = -(log(0.5)+log(1.0-exp(-2.0*rargs.indel_rate*rargs.branch_length))
+ + log(rargs.avgaln)-log(rargs.avgaln+1.0)
+ - log(rargs.indel_mean-1.0));
+ dB = -(log(rargs.indel_mean-1.0)-log(rargs.indel_mean));
+ dC = 0.0;
+ dF = dH = 0.0;
+ dG = 0.0;
+
+ return true;
+}
+
+
+
diff --git a/src/models.h b/src/models.h
new file mode 100644
index 0000000..7130113
--- /dev/null
+++ b/src/models.h
@@ -0,0 +1,116 @@
+/****************************************************************************
+ * Copyright (C) 2005-2010 Reed A. Cartwright, PhD <reed at scit.us> *
+ * *
+ * This program is free software: you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation, either version 3 of the License, or *
+ * (at your option) any later version. *
+ * *
+ * This program is distributed in the hope that it will be useful, *
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of *
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
+ * GNU General Public License for more details. *
+ * *
+ * You should have received a copy of the GNU General Public License *
+ * along with this program. If not, see <http://www.gnu.org/licenses/>. *
+ ****************************************************************************/
+
+#ifndef MODELS_H
+#define MODELS_H
+
+#include <cmath>
+#include <cfloat>
+#include <limits>
+
+#include <boost/cstdint.hpp>
+#include <boost/preprocessor/cat.hpp>
+
+#include "ngila_app.h"
+#include "matparser.h"
+#include "seqdb.h"
+
+#ifndef TABLE_CELL_BITSIZE
+# define TABLE_CELL_BITSIZE 16
+#endif
+
+typedef BOOST_PP_CAT(boost::int, BOOST_PP_CAT(TABLE_CELL_BITSIZE, _t)) _travel_cell;
+
+struct cost_model
+{
+ double dA, dB, dC;
+ double dF, dG, dH;
+ sub_matrix mCost;
+
+ virtual bool create(const ngila_app::args &rargs);
+ virtual double offset(const std::string &seqA,
+ const std::string &seqB) const
+ {
+ return 0.0;
+ }
+
+ inline double gapcost(size_t L)
+ {
+ if(L < 1)
+ return 0.0;
+ if(L > (size_t)std::numeric_limits<_travel_cell>::max())
+ return std::numeric_limits<double>::max()/16.0;
+ return dA+dB*L+dC*log((double)L);
+ }
+ inline double freegapcost(size_t L)
+ {
+ if(L < 1)
+ return 0.0;
+ if(L > (size_t)std::numeric_limits<_travel_cell>::max())
+ return std::numeric_limits<double>::max()/16.0;
+ return dF+dG*L+dH*log((double)L);
+ }
+ inline size_t kstar(double x, double y)
+ {
+ return static_cast<size_t>(floor(x/(1.0-exp((-y+dB*x)/dC)))); //needs to handle dC == 0.0?
+ }
+ inline size_t kstar_f(double x, double y)
+ {
+ return static_cast<size_t>(floor(x/(1.0-exp((-y+dG*x)/dH)))); //needs to handle dC == 0.0?
+ }
+};
+
+struct k2p_model : public cost_model {
+ double dEnd;
+
+ virtual bool create(const ngila_app::args &rargs);
+ virtual double offset(const std::string &seqA,
+ const std::string &seqB) const;
+protected:
+ double dNucScale;
+};
+
+struct zeta_model : public k2p_model {
+ virtual bool create(const ngila_app::args &rargs);
+};
+
+struct geo_model : public k2p_model {
+ virtual bool create(const ngila_app::args &rargs);
+};
+
+struct aa_model : public cost_model {
+ double dEnd;
+ double vAACost[128];
+
+ virtual bool create(const ngila_app::args &rargs);
+ virtual double offset(const std::string &seqA,
+ const std::string &seqB) const;
+protected:
+ double dAaScale;
+};
+
+struct aazeta_model : public aa_model {
+ virtual bool create(const ngila_app::args &rargs);
+};
+
+struct aageo_model : public aa_model {
+ virtual bool create(const ngila_app::args &rargs);
+};
+
+
+#endif
+
diff --git a/src/ngila-256.png b/src/ngila-256.png
new file mode 100644
index 0000000..337f679
Binary files /dev/null and b/src/ngila-256.png differ
diff --git a/src/ngila-256.psd b/src/ngila-256.psd
new file mode 100644
index 0000000..1fd3f31
Binary files /dev/null and b/src/ngila-256.psd differ
diff --git a/src/ngila.cmds b/src/ngila.cmds
new file mode 100644
index 0000000..f50ed77
--- /dev/null
+++ b/src/ngila.cmds
@@ -0,0 +1,69 @@
+/****************************************************************************
+ * Copyright (C) 2005-2010 Reed A. Cartwright, PhD <reed at scit.us> *
+ * *
+ * This program is free software: you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation, either version 3 of the License, or *
+ * (at your option) any later version. *
+ * *
+ * This program is distributed in the hope that it will be useful, *
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of *
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
+ * GNU General Public License for more details. *
+ * *
+ * You should have received a copy of the GNU General Public License *
+ * along with this program. If not, see <http://www.gnu.org/licenses/>. *
+ ****************************************************************************/
+
+#include "xm.h"
+
+/***************************************************************************
+ * X-Macro List *
+ * *
+ * XCMD(lname, sname, desc, type, def) *
+ ***************************************************************************/
+
+XCMD((model), (m), "alignment model: zeta|geo|aazeta|aageo|cost", std::string, std::string("zeta"))
+XCMD((branch)(length), (t), "sequence separation time", double, _DL(0.1, "0.1"))
+XCMD((ratio), (k), "transition-transversion ratio", double, 2.0)
+XCMD((indel)(rate), (r), "relative rate of insertion and deletion", double, _DL(0.1, "0.1"))
+XCMD((indel)(slope), (z), "slope param of insertion and deletion (for zeta)", double, 1.7)
+XCMD((indel)(mean), (q), "mean of insertion and deletion (for geo)", double, 10.0)
+XCMD((avgaln), (l), "average alignment skeleton length", double, 1000)
+
+XCMD((cost)(match), (i), "residue match cost", double, 0.0)
+XCMD((cost)(mismatch), (j), "residue mismatch cost", double, 1.0)
+XCMD((cost)(matrix), (x), "read residue cost matrix from file", std::string, std::string(""))
+XCMD((cost)(intersection), (a), "gap intersection cost", double, 1.0)
+XCMD((cost)(linear), (b), "gap linear cost", double, 1.0)
+XCMD((cost)(logarithmic), (c), "gap logarithmic cost", double, 0.0)
+XCMD((cost)(intersection)(free), (f), "gap intersection cost", double, 0.0)
+XCMD((cost)(linear)(free), (g), "gap linear cost", double, 0.0)
+XCMD((cost)(logarithmic)(free), (h), "gap logarithmic cost", double, 0.0)
+
+XCMD((free)(end)(gaps), (e), "free end gaps", bool, false)
+XCMD((case)(insensitivity), (I), "sequences are case insensitive", bool, false)
+XCMD((threshold)(larger), (M), "threshold for O(MN) alignment", size_t, 10000)
+XCMD((threshold)(smaller), (N), "threshold for O(MN) alignment", size_t, 10000)
+XCMD((remove)(gaps), (G), "realign after removing specific gap characters", std::string, std::string("-+="))
+#ifdef USE_THREADS
+XCMD((threads), (T), "number of worker threads to use", size_t, 1)
+#endif
+XCMD((output), (o), "output file", std::string, std::string("-"))
+
+XCMD((pairs), , "control how pairs are determined: first|all|each", std::string, std::string("first"))
+XCMD((no)(scale), , "do not force identical sequences to have a cost of zero", bool, false)
+XCMD((no)(header), , "exclude a header of sequence names in table output format", bool, false)
+XCMD((const)(align), , "sort and align sequences regardless of direction and order", size_t, 0)
+XCMD((separator), , "text to use as a separator of alignments", std::string, std::string(""))
+XCMD((quiet), , "disable all warnings and error messages", bool, false)
+XCMD((ngilarc), , "alternative path for rc file", std::string, std::string(""))
+XCMD((version), , "display version information", bool, false)
+XCMD((help), , "display help message", bool, false)
+XCMD((desktop), , "behave as if called on the desktop", bool, false)
+
+
+/***************************************************************************
+ * cleanup *
+ ***************************************************************************/
+#include "xm.h"
diff --git a/src/ngila.cpp b/src/ngila.cpp
new file mode 100644
index 0000000..250fd94
--- /dev/null
+++ b/src/ngila.cpp
@@ -0,0 +1,437 @@
+/****************************************************************************
+ * Copyright (C) 2005-2010 Reed A. Cartwright, PhD <reed at scit.us> *
+ * *
+ * This program is free software: you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation, either version 3 of the License, or *
+ * (at your option) any later version. *
+ * *
+ * This program is distributed in the hope that it will be useful, *
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of *
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
+ * GNU General Public License for more details. *
+ * *
+ * You should have received a copy of the GNU General Public License *
+ * along with this program. If not, see <http://www.gnu.org/licenses/>. *
+ ****************************************************************************/
+
+#include "ngila.h"
+#include <sstream>
+#include <iomanip>
+#include <iostream>
+#include <limits>
+
+#include <boost/preprocessor.hpp>
+#include <boost/foreach.hpp>
+#include <boost/config.hpp>
+#include <boost/scoped_ptr.hpp>
+//BOOST_DISABLE_ASSERTS
+#include <boost/multi_array.hpp>
+#include <boost/algorithm/string/predicate.hpp>
+
+#ifdef USE_THREADS
+# include <boost/thread.hpp>
+#endif
+
+#include "ngila_app.h"
+#include "seqdb.h"
+#include "matparser.h"
+#include "models.h"
+#include "align.h"
+
+#define foreach BOOST_FOREACH
+
+using namespace std;
+
+typedef std::pair<size_t,size_t> job_pair;
+typedef std::vector<job_pair> pair_vec;
+
+#ifdef USE_THREADS
+struct job_info {
+
+ typedef boost::mutex::scoped_lock lock;
+
+ const seq_db &db;
+ const pair_vec &pairs;
+ const cost_model &model;
+ const size_t const_align;
+ const int out_format;
+ const int direction;
+
+ job_info(const seq_db &d, const pair_vec &v, const cost_model &m, size_t ca,
+ int of, int dir, ostream &o)
+ : db(d), pairs(v), model(m), const_align(ca), out_format(of),
+ direction((ca&16)?0:dir),
+ pos(0), out(o) { }
+
+ size_t next_pair(job_pair &p) {
+ lock lk(np_monitor);
+ if(pos >= pairs.size())
+ return size_t(-1);
+ p = pairs[pos];
+ return pos++;
+ }
+
+ void print_aln(const alignment &aln, double dcost, double dident, bool swapped) {
+ lock lk(print_monitor);
+ // round to nearest billionth
+ dcost = floor(dcost*1e9+0.5)/1e9;
+ dident = floor(dident*1e9+0.5)/1e9;
+
+ if(out_format < 3)
+ aln.print(out, out_format, dcost, direction, swapped);
+
+ }
+
+
+private:
+ pair_vec::size_type pos;
+ boost::mutex np_monitor, print_monitor;
+ ostream &out;
+};
+
+
+struct ngila_worker {
+ ngila_worker(const aligner &a, job_info &i) :
+ alner(a), info(i) {
+ }
+ void operator()() {
+ // swap a and b so that a's hashed position is lower
+ job_pair p;
+ size_t pos;
+ while((pos = info.next_pair(p)) != size_t(-1)) {
+ size_t a = p.first;
+ size_t b = p.second;
+ bool swapped = false;
+ if(!(info.const_align & 8) && info.db.db().project<hashid>(info.db.db().begin()+a) >
+ info.db.db().project<hashid>(info.db.db().begin()+b) ) {
+ swap(a,b);
+ swapped = true;
+ }
+
+ alignment aln(info.db[a], info.db[b]);
+ double dcost = alner.align(aln);
+ dcost += info.model.offset(info.db[a].dna, info.db[b].dna);
+ double dident = aln.identity();
+
+ info.print_aln(aln, dcost, dident, swapped);
+ }
+ }
+private:
+ aligner alner;
+ job_info &info;
+};
+#endif
+
+
+int main(int argc, char *argv[])
+{
+ int ret = EXIT_FAILURE;
+ try {
+ ngila_app app(argc, argv);
+ ret = app.run();
+ } catch(exception &e) {
+ CERROR(e.what());
+ }
+ return ret;
+}
+
+ngila_app::ngila_app(int argc, char* argv[]) : desc("Allowed Options") {
+ runname = argv[0];
+ try {
+ desc.add_options()
+ #define XCMD(lname, sname, desc, type, def) ( \
+ _S(lname) _IFD(sname, "," BOOST_PP_STRINGIZE sname), \
+ po::value< type >(&arg._V(lname))->default_value(def), \
+ desc )
+ #include "ngila.cmds"
+ #undef XCMD
+ ;
+ indesc.add_options()("input", po::value< vector<string> >(&arg.input), "input files");
+ indesc.add(desc);
+ pdesc.add("input", -1);
+ po::store(po::command_line_parser(argc, argv).options(indesc).positional(pdesc).run(), vm);
+ po::notify(vm);
+ // check to see if ngilarc is specified as stdin
+ if(arg.ngilarc == "-") {
+ po::store(po::parse_config_file(cin, desc), vm);
+ po::notify(vm);
+ } else {
+ std::ifstream ifs;
+ // if ngilarc is specified try to open it
+ if(!arg.ngilarc.empty()) {
+ ifs.open(arg.ngilarc.c_str());
+ if(!ifs.is_open()) {
+ string sse = "unable to open argument file ";
+ sse += arg.ngilarc;
+ throw std::runtime_error(sse);
+ }
+ } else {
+ // if ngilarc is not specified, check for it in HOME
+#ifdef BOOST_WINDOWS
+ char buf[_MAX_PATH];
+ _makepath_s(buf, getenv("HOMEDRIVE"), getenv("HOMEPATH"), "ngilarc", "txt");
+ ifs.open(buf);
+#else
+ char *h;
+ if((h = getenv("HOME")) != NULL) {
+ string buf = h;
+ buf += "/.ngilarc";
+ ifs.open(buf.c_str());
+ }
+#endif
+ }
+ if(ifs.is_open()) {
+ po::store(po::parse_config_file(ifs, desc), vm);
+ po::notify(vm);
+ }
+ }
+ } catch (exception &e) {
+ CERROR(e.what());
+ throw std::runtime_error("unable to process command line");
+ }
+}
+
+template<std::size_t _N>
+std::size_t key_switch(const std::string &ss, const std::string (&key)[_N]) {
+ using boost::algorithm::starts_with;
+ for(std::size_t i=0;i<_N;++i) {
+ if(starts_with(key[i], ss))
+ return i;
+ }
+ return (std::size_t)-1;
+}
+
+int ngila_app::run()
+{
+ if(arg.version) {
+ cerr << endl << VERSION_MSG << endl << endl;
+ return EXIT_SUCCESS;
+ }
+ if(arg.help || arg.input.empty()) {
+ cerr << endl << VERSION_MSG << endl << endl;
+ cerr << "Usage:\n "
+ << runname << " [options] input.fas ..."
+ << endl << endl;
+ cerr << desc << endl;
+ return EXIT_SUCCESS;
+ }
+ if(arg.quiet)
+ cerr.clear(ios::failbit);
+
+ seq_db mydb(arg.remove_gaps);
+ for(vector<string>::const_iterator cit = arg.input.begin(); cit != arg.input.end(); ++cit) {
+ if(!mydb.parse_file(cit->c_str(), true, arg.case_insensitivity)) {
+ CERROR("parsing of \'" << cit->c_str() << "\' failed.");
+ return EXIT_FAILURE;
+ }
+ }
+
+ int direction = 0;
+ if(arg.const_align != 0) {
+ direction = mydb.unique_sort(seq_db::DIR_ORI | (arg.const_align & 7));
+ mydb.transform(direction);
+ if(!(arg.const_align & 8)) { // return to original order
+ mydb.rearrange(mydb.db().get<id>().begin());
+ }
+ }
+
+
+ boost::scoped_ptr<cost_model> pmod;
+ string model_keys[] = {
+ string("zeta"), string("geo"),
+ string("aazeta"), string("aageo"),
+ string("cost")
+ };
+ switch(key_switch(arg.model, model_keys))
+ {
+ case 0:
+ pmod.reset(new zeta_model);
+ break;
+ case 1:
+ pmod.reset(new geo_model);
+ break;
+ case 2:
+ pmod.reset(new aazeta_model);
+ break;
+ case 3:
+ pmod.reset(new aageo_model);
+ break;
+ case 4:
+ pmod.reset(new cost_model);
+ break;
+ default:
+ CERROR("unknown model \'" << arg.model << "\'.");
+ return EXIT_FAILURE;
+ };
+
+ if(!pmod->create(arg))
+ {
+ CERROR("creating model \'" << arg.model << "\'.");
+ return EXIT_FAILURE;
+ }
+ if(mydb.size() < 2)
+ {
+ CERROR("two or more sequences are required for alignment.");
+ return EXIT_FAILURE;
+ }
+ pair_vec pvec;
+ string pairs_keys[] = { string("first"), string("all"), string("each") };
+ seq_db::size_type table_size = mydb.size();
+
+ switch(key_switch(arg.pairs, pairs_keys)) {
+ case 0:
+ pvec.push_back(make_pair(0,1));
+ table_size = 2;
+ break;
+ case 1:
+ for(size_t i=0;i<mydb.size();++i)
+ for(size_t j=i+1;j<mydb.size();++j)
+ pvec.push_back(make_pair(i,j));
+ break;
+ case 2:
+ for(size_t i=0;i<mydb.size()-1;i+=2)
+ pvec.push_back(make_pair(i,i+1));
+ break;
+ default:
+ CERROR("unknown pairs option \'" << arg.pairs << "\'.");
+ return EXIT_FAILURE;
+ };
+
+ string format_keys[] = {
+ string("aln"), string("fasta"), string("phylip"),
+ string("dist"), string("dist-c"),
+ string("dist-d"), string("dist-i")
+ };
+ int out_format = 0;
+ if(!arg.output.empty()) {
+ string ssFormat;
+ string::size_type pos = arg.output.find_first_of(':');
+#ifdef BOOST_WINDOWS
+ if(pos != string::npos && pos != 1) {
+#else
+ if(pos != string::npos) { // format:file
+#endif
+ ssFormat = arg.output.substr(0, pos);
+ arg.output.erase(0, pos+1);
+ } else { // file.format
+ pos = arg.output.find_last_of('.');
+ if(pos != string::npos)
+ ssFormat = arg.output.substr(pos+1);
+ }
+ if(!ssFormat.empty()) {
+ out_format = key_switch(ssFormat, format_keys);
+ if(out_format == -1) {
+ CERROR("unknown output format \'" << ssFormat << "\'.");
+ return EXIT_FAILURE;
+ }
+ }
+ }
+
+ ofstream fout;
+ if(!(arg.output.empty() || arg.output == "-")) {
+ fout.open(arg.output.c_str(), ios_base::out|ios_base::trunc);
+ if(!fout.is_open()) {
+ CERROR("unable to open output file \'" << arg.output << "\'.");
+ return EXIT_FAILURE;
+ }
+ }
+ ostream &myout = fout.is_open() ? static_cast<ostream&>(fout)
+ : static_cast<ostream&>(cout);
+ typedef boost::multi_array<float, 2> dist_mat;
+ dist_mat dist_table;
+ const bool do_dist = (out_format >= 3);
+ if(do_dist) {
+ dist_table.resize(boost::extents[table_size][table_size]);
+ fill(dist_table.data(), dist_table.data()+dist_table.num_elements(), numeric_limits<float>::quiet_NaN());
+ float diag = (out_format == 6) ? 1.0f : 0.0f;
+ for(seq_db::size_type i = 0; i < table_size; ++i)
+ dist_table[i][i] = diag;
+ }
+
+ aligner alner(*pmod, arg.threshold_larger, arg.threshold_smaller, arg.free_end_gaps);
+
+#ifdef USE_THREADS
+ boost::thread_group threads;
+ job_info info(mydb, pvec, *pmod, arg.const_align,
+ out_format, direction, myout);
+
+ for (size_t i = 0; i < arg.threads; ++i)
+ threads.create_thread(ngila_worker(alner,info));
+ threads.join_all();
+
+#else
+
+ for(pair_vec::const_iterator cit = pvec.begin(); cit != pvec.end(); ++cit) {
+ if(!do_dist && cit != pvec.begin())
+ myout << arg.separator << endl;
+ // swap a and b so that a's hashed position is lower
+ size_t a = cit->first;
+ size_t b = cit->second;
+ bool swapped = false;
+ if(!(arg.const_align & 8) && mydb.db().project<hashid>(mydb.db().begin()+a) >
+ mydb.db().project<hashid>(mydb.db().begin()+b) ) {
+ swap(a,b);
+ swapped = true;
+ }
+
+ alignment aln(mydb[a], mydb[b]);
+ double dcost = alner.align(aln);
+ dcost += pmod->offset(mydb[a].dna, mydb[b].dna);
+ double dident = aln.identity();
+
+ // round to nearest billionth
+ dcost = floor(dcost*1e9+0.5)/1e9;
+ dident = floor(dident*1e9+0.5)/1e9;
+
+ if(!do_dist) {
+ aln.print(myout, out_format, dcost,
+ ((arg.const_align & 16) ? 0 : direction), swapped);
+ } else {
+ double upper,lower;
+ if(a > b)
+ swap(a,b);
+ switch(out_format) {
+ default:
+ case 3: //dist
+ lower = 1.0-dident;
+ upper = dcost;
+ break;
+ case 4: //dist-c
+ upper = lower = dcost;
+ break;
+ case 5: //dist-d
+ upper = lower = 1.0-dident;
+ break;
+ case 6: //dist-i
+ upper = lower = dident;
+ break;
+ }
+ dist_table[a][b] = (float)upper;
+ dist_table[b][a] = (float)lower;
+ }
+ }
+ if(do_dist) {
+ myout << setprecision(10);
+ if(!arg.no_header) {
+ for(seq_db::size_type i = 0; i < table_size-1; ++i)
+ myout << mydb[i].name << "\t";
+ myout << mydb[table_size-1].name << endl;
+ }
+ for(seq_db::size_type i = 0; i < table_size; ++i) {
+ for(seq_db::size_type j = 0; j < table_size-1; ++j)
+ myout << dist_table[i][j] << "\t";
+ myout << dist_table[i][table_size-1] << endl;
+ }
+ }
+
+#endif
+
+ if(arg.desktop) {
+ cerr << "\nPress any key to continue." << endl;
+ char x;
+ cin.get(x);
+ }
+
+ return EXIT_SUCCESS;
+}
diff --git a/src/ngila.desktop b/src/ngila.desktop
new file mode 100644
index 0000000..08da7e8
--- /dev/null
+++ b/src/ngila.desktop
@@ -0,0 +1,13 @@
+[Desktop Entry]
+Version=1.0
+Encoding=UTF-8
+Name=Ngila
+Comment=Log-Affine Alignments
+Exec=ngila --desktop %f
+TryExec=ngila
+Icon=ngila
+Terminal=true
+X-MultipleArgs=false
+Type=Application
+Categories=Application;Education;Science;Biology;
+StartupNotify=false
diff --git a/src/ngila.h b/src/ngila.h
new file mode 100644
index 0000000..7217fa3
--- /dev/null
+++ b/src/ngila.h
@@ -0,0 +1,73 @@
+/****************************************************************************
+ * Copyright (C) 2005-2010 Reed A. Cartwright, PhD <reed at scit.us> *
+ * *
+ * This program is free software: you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation, either version 3 of the License, or *
+ * (at your option) any later version. *
+ * *
+ * This program is distributed in the hope that it will be useful, *
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of *
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
+ * GNU General Public License for more details. *
+ * *
+ * You should have received a copy of the GNU General Public License *
+ * along with this program. If not, see <http://www.gnu.org/licenses/>. *
+ ****************************************************************************/
+
+#ifndef NGILA_H
+#define NGILA_H
+
+#ifdef HAVE_CONFIG_H
+# include "config.h"
+#endif
+
+#include <iostream>
+#include <vector>
+#include <deque>
+#include <map>
+#include <string>
+#include <cstring>
+#include <fstream>
+#include <ostream>
+#include <iterator>
+
+#define VERSION_MSG PACKAGE_STRING "\n" \
+ " Copyright (C) 2005-2010 Reed A. Cartwright, PhD <reed at scit.us>\n" \
+ " Report bugs to " PACKAGE_BUGREPORT
+
+// Utility functions
+#define CERRORR(err_msg) ((std::cerr << "ERROR: " << err_msg << std::endl), false);
+#define CERROR(err_msg) (std::cerr << "ERROR: " << err_msg << std::endl);
+
+namespace std {
+template<typename _Tp, typename _CharT, typename _Traits>
+basic_ostream<_CharT, _Traits>&
+operator<<(basic_ostream<_CharT, _Traits>& os, const std::vector<_Tp> &v)
+{
+ if(v.size() == 1)
+ os << v.front();
+ else if(v.size() > 1)
+ {
+ std::copy(v.begin(), v.end()-1, std::ostream_iterator<_Tp, _CharT, _Traits>(os, " "));
+ os << v.back();
+ }
+ return os;
+}
+template<typename _Tp, typename _CharT, typename _Traits>
+basic_ostream<_CharT, _Traits>&
+operator<<(basic_ostream<_CharT, _Traits>& os, const std::deque<_Tp> &v)
+{
+ if(v.size() == 1)
+ os << v.front();
+ else if(v.size() > 1)
+ {
+ std::copy(v.begin(), v.end()-1, std::ostream_iterator<_Tp, _CharT, _Traits>(os, " "));
+ os << v.back();
+ }
+ return os;
+}
+}
+
+#endif //NGILA_H
+
diff --git a/src/ngila.ico b/src/ngila.ico
new file mode 100644
index 0000000..543568e
Binary files /dev/null and b/src/ngila.ico differ
diff --git a/src/ngila.png b/src/ngila.png
new file mode 100644
index 0000000..401666c
Binary files /dev/null and b/src/ngila.png differ
diff --git a/src/ngila.rc b/src/ngila.rc
new file mode 100644
index 0000000..df58e49
--- /dev/null
+++ b/src/ngila.rc
@@ -0,0 +1 @@
+IDI_ICON1 ICON DISCARDABLE "ngila.ico"
diff --git a/src/ngila_app.h b/src/ngila_app.h
new file mode 100644
index 0000000..270fd99
--- /dev/null
+++ b/src/ngila_app.h
@@ -0,0 +1,111 @@
+/****************************************************************************
+ * Copyright (C) 2005-2010 Reed A. Cartwright, PhD <reed at scit.us> *
+ * *
+ * This program is free software: you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation, either version 3 of the License, or *
+ * (at your option) any later version. *
+ * *
+ * This program is distributed in the hope that it will be useful, *
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of *
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
+ * GNU General Public License for more details. *
+ * *
+ * You should have received a copy of the GNU General Public License *
+ * along with this program. If not, see <http://www.gnu.org/licenses/>. *
+ ****************************************************************************/
+
+#ifndef NGILA_APP_H
+#define NGILA_APP_H
+#pragma once
+
+#include "ngila.h"
+
+#include <vector>
+#include <utility>
+
+#include <boost/logic/tribool.hpp>
+#include <boost/logic/tribool_io.hpp>
+#include <boost/program_options.hpp>
+namespace po = boost::program_options;
+
+/****************************************************************************
+ * class ngila_app *
+ ****************************************************************************/
+
+class ngila_app {
+public:
+ ngila_app(int argc, char *argv[]);
+ virtual ~ngila_app() { }
+ virtual int run();
+
+ std::string runname;
+ po::options_description desc, indesc;
+ po::positional_options_description pdesc;
+ po::variables_map vm;
+
+ struct args {
+ // use X-Macros to specify argument variables
+ #define XCMD(lname, sname, desc, type, def) type _V(lname) ;
+ #include "ngila.cmds"
+ #undef XCMD
+ std::vector< std::string > input;
+ };
+ enum { CONST_ALIGN_ORDER = 8, CONST_ALIGN_SEQS = 16 };
+
+protected:
+ args arg;
+
+};
+
+namespace boost { namespace program_options {
+template<>
+inline typed_value<bool>* value(bool* v) {
+ return bool_switch(v);
+}
+template<>
+inline typed_value<boost::tribool>* value(boost::tribool* v) {
+ typed_value<boost::tribool>* r = new typed_value<boost::tribool>(v);
+ r->implicit_value(true, "on");
+ return r;
+}
+
+// modified from boost/libs/program_options/src/value_semantic.cpp
+inline void validate(boost::any& v, const std::vector<std::string>& xs, boost::tribool*, int) {
+ check_first_occurrence(v);
+ std::string s(get_single_string(xs, true));
+
+ for (size_t i = 0; i < s.size(); ++i)
+ s[i] = char(tolower(s[i]));
+
+ if (s.empty() || s == "on" || s == "yes" || s == "1" || s == "true")
+ v = any(boost::tribool(true));
+ else if (s == "off" || s == "no" || s == "0" || s == "false")
+ v = any(boost::tribool(false));
+ else if (s == "null" || s == "maybe" || s == "2" || s == "indeterminate")
+ v = any(boost::tribool(boost::indeterminate));
+ else
+ boost::throw_exception(invalid_option_value(s));
+}
+#if !defined(BOOST_NO_STD_WSTRING)
+inline void validate(any& v, const std::vector<std::wstring>& xs, boost::tribool*, int)
+{
+ check_first_occurrence(v);
+ std::wstring s(get_single_string(xs, true));
+
+ for (size_t i = 0; i < s.size(); ++i)
+ s[i] = wchar_t(tolower(s[i]));
+
+ if (s.empty() || s == L"on" || s == L"yes" || s == L"1" || s == L"true")
+ v = any(boost::tribool(true));
+ else if (s == L"off" || s == L"no" || s == L"0" || s == L"false")
+ v = any(boost::tribool(false));
+ else if (s == L"null" || s == L"maybe" || s == L"2" || s == L"indeterminate")
+ v = any(boost::tribool(boost::indeterminate));
+ else
+ boost::throw_exception(invalid_option_value(s));
+}
+#endif
+}}
+
+#endif
diff --git a/src/seqdb.cpp b/src/seqdb.cpp
new file mode 100644
index 0000000..087a80e
--- /dev/null
+++ b/src/seqdb.cpp
@@ -0,0 +1,147 @@
+/****************************************************************************
+ * Copyright (C) 2005-2010 Reed A. Cartwright, PhD <reed at scit.us> *
+ * *
+ * This program is free software: you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation, either version 3 of the License, or *
+ * (at your option) any later version. *
+ * *
+ * This program is distributed in the hope that it will be useful, *
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of *
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
+ * GNU General Public License for more details. *
+ * *
+ * You should have received a copy of the GNU General Public License *
+ * along with this program. If not, see <http://www.gnu.org/licenses/>. *
+ ****************************************************************************/
+
+#include "ngila.h"
+
+#include <cstring>
+#include <string>
+#include <algorithm>
+#include <boost/bind.hpp>
+#include <boost/foreach.hpp>
+#include <boost/functional/hash.hpp>
+
+#include "seqdb.h"
+#include "seqparser.h"
+
+#define foreach BOOST_FOREACH
+
+using namespace std;
+using namespace boost;
+using namespace boost::multi_index;
+
+
+bool seq_db::parse_file(const char *csfile, bool bappend, bool bi)
+{
+ if(!bappend)
+ clear();
+
+ stack<string> my_stack;
+ size_type pos;
+ seq_grammar my_grammar(my_stack, *this,pos);
+
+ if(strcmp(csfile, "-")==0) {
+ string ss;
+ getline(cin, ss, '\004'); // EOT, end of transmission, ^D
+ if(ss.empty())
+ return CERRORR("unable to open stdin.");
+
+ parse_info<string::const_iterator> info = parse(ss.begin(), ss.end(), my_grammar);
+ if (!info.full)
+ return CERRORR("unable to parse stdin.");
+ } else {
+ file_iterator<char> file_first(csfile);
+ if(!file_first)
+ return CERRORR("unable to open \'" << csfile << "\'");
+ file_iterator<char> file_last = file_first.make_end();
+
+ parse_info< file_iterator<char> > info = parse(file_first, file_last, my_grammar);
+ if (!info.full)
+ return CERRORR("unable to parse \'" << csfile << "\'");
+ }
+
+ for(container::iterator it = cont.begin(); it != cont.end(); ++it)
+ cont.modify(it, boost::bind(&seq_data::sanitize, _1, bi));
+
+ return true;
+}
+
+struct hash_data {
+ std::size_t hashed;
+ int dir;
+ typedef boost::reference_wrapper<const seq_data> wrap;
+ wrap ref;
+
+ hash_data(std::size_t h, int d, const seq_data &s) :
+ hashed(h), dir(d), ref(s) {
+ }
+ operator const seq_data&() const {
+ return ref;
+ }
+};
+
+struct full_key : composite_key< hash_data,
+ member< hash_data, int, &hash_data::dir>,
+ member< hash_data, size_t, &hash_data::hashed>
+>{};
+
+typedef multi_index_container< hash_data, indexed_by<
+ ordered_non_unique<
+ member< hash_data, size_t, &hash_data::hashed>, greater<size_t> >,
+ ordered_non_unique<
+ full_key, composite_key_result_greater<full_key::result_type> >
+ > > hash_container;
+
+typedef hash_container::nth_index<1>::type hc_by_both;
+
+int seq_db::unique_sort(int directions) {
+ seq_data seq("a", "aa");
+ hash_data hd(0, 0, seq);
+ const seq_data &ww = (hd);
+
+ hash_container hc;
+ foreach(const seq_data & sd, cont) {
+ if(directions & DIR_ORI) // hash sequence
+ hc.insert(hash_data(hash_range(sd.dna.begin(), sd.dna.end()), DIR_ORI, sd));
+ if(directions & DIR_REV) // hash reverse sequence
+ hc.insert(hash_data(hash_range(sd.dna.rbegin(), sd.dna.rend()), DIR_REV, sd));
+ if(directions & DIR_COM) // hash complement sequence
+ hc.insert(hash_data(hash_range(complementer(sd.dna.begin()),
+ complementer(sd.dna.end())), DIR_COM, sd));
+ if(directions & DIR_RVC) // hash reverse complement
+ hc.insert(hash_data(hash_range(complementer(sd.dna.rbegin()),
+ complementer(sd.dna.rend())), DIR_RVC, sd));
+ }
+ size_t hh = hc.begin()->hashed;
+ int dd = hc.begin()->dir, mask = -1;
+ foreach(const hash_data &h, hc) {
+ if(h.hashed == hh) {
+ dd |= h.dir;
+ continue;
+ }
+ dd &= mask;
+ switch(dd) {
+ case DIR_ORI: case DIR_REV: case DIR_COM: case DIR_RVC:
+ break;
+ default:
+ mask = dd;
+ case 0:
+ hh = h.hashed;
+ dd = h.dir;
+ continue;
+ }
+ break;
+ }
+ // dd now represents the best direction(s)
+ // use an array to break any ties 8 > 1 > 2 > 4
+ static const int ties[] = { 8,1,2,1,4,1,2,1,8,8,8,8,8,8,8,8 };
+ dd = ties[dd&15];
+ // sort
+ cont.rearrange(hc.get<1>().find(dd));
+ // store this order
+ cont.get<hashid>().rearrange(cont.begin());
+ return dd;
+}
diff --git a/src/seqdb.h b/src/seqdb.h
new file mode 100644
index 0000000..8bca23f
--- /dev/null
+++ b/src/seqdb.h
@@ -0,0 +1,190 @@
+/****************************************************************************
+ * Copyright (C) 2005-2010 Reed A. Cartwright, PhD <reed at scit.us> *
+ * *
+ * This program is free software: you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation, either version 3 of the License, or *
+ * (at your option) any later version. *
+ * *
+ * This program is distributed in the hope that it will be useful, *
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of *
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
+ * GNU General Public License for more details. *
+ * *
+ * You should have received a copy of the GNU General Public License *
+ * along with this program. If not, see <http://www.gnu.org/licenses/>. *
+ ****************************************************************************/
+
+#ifndef SEQDB_H
+#define SEQDB_H
+
+#include <string>
+#include <vector>
+#include <map>
+#include <algorithm>
+
+#include <boost/bind.hpp>
+#include <boost/lambda/lambda.hpp>
+#include <boost/multi_index_container.hpp>
+#include <boost/multi_index/ordered_index.hpp>
+#include <boost/multi_index/sequenced_index.hpp>
+#include <boost/multi_index/random_access_index.hpp>
+#include <boost/multi_index/identity.hpp>
+#include <boost/multi_index/composite_key.hpp>
+#include <boost/multi_index/member.hpp>
+#include <boost/functional/hash.hpp>
+
+
+inline char complement(char x) {
+ static char dnac[] = "\0\x01\x02\x03\x04\x05\x06\a\b\t\n\x0b\f\r\x0e\x0f"
+ "\x10\x11\x12\x13\x14\x15\x16\x17\x18\x19\x1a\x1b\x1c\x1d\x1e\x1f"
+ " !\"#$%&'()*+,-./0123456789:;<=>?@"
+ "TVGHEFCDIJMLKNOPQYSAABWXRZ" "[\\]^_`"
+ "tvghefcdijmlknopqysaabwxrz" "{|}~\x7f";
+
+ return dnac[x];
+}
+
+template<class _Iterator>
+class _complementer : public _Iterator {
+public:
+ typedef _Iterator base_type;
+ typedef typename base_type::value_type value_type;
+
+ value_type operator*() const {
+ return value_type(complement(**static_cast<const base_type*>(this)));
+ }
+
+ _complementer(base_type b) : base_type(b) { }
+};
+
+template<class _Iterator>
+_complementer<_Iterator> complementer(_Iterator it) {
+ return _complementer<_Iterator>(it);
+}
+
+class seq_db;
+
+struct seq_data {
+ std::string name;
+ std::string dna;
+
+ seq_data(const std::string &n, const std::string &s) : name(n), dna(s) {
+ }
+
+ // append to dna
+ inline seq_data& operator+=(const seq_data& a) {
+ dna += a.dna;
+ return *this;
+ }
+
+ // remove characters in dna that belong to g
+ inline seq_data& operator-=(const std::string& g) {
+ std::string::size_type (std::string::*gf)(std::string::value_type, std::string::size_type) const
+ = &std::string::find;
+ dna.erase( remove_if(dna.begin(), dna.end(),
+ boost::bind(gf, g, _1, 0) != std::string::npos),
+ dna.end());
+ return *this;
+ }
+
+ inline void sanitize(bool bi=false) {
+ std::replace_if(name.begin(), name.end(), std::ptr_fun(::isspace), '_');
+ if(bi)
+ std::transform(dna.begin(), dna.end(), dna.begin(), std::ptr_fun(::toupper));
+ }
+
+ inline void transform(int dir);
+};
+
+//tags
+struct name {};
+struct id {};
+struct hashid {};
+
+class seq_db {
+public:
+ typedef boost::multi_index_container< seq_data, boost::multi_index::indexed_by<
+ boost::multi_index::random_access<>, // current order
+ boost::multi_index::random_access<boost::multi_index::tag<id> >, // original order
+ boost::multi_index::random_access<boost::multi_index::tag<hashid> >, // hashed order
+ boost::multi_index::ordered_unique<
+ boost::multi_index::tag<name>,
+ boost::multi_index::member< seq_data, std::string, &seq_data::name> >
+ > > container;
+ typedef container::size_type size_type;
+
+ enum {
+ DIR_REV = 1, DIR_COM = 2, DIR_RVC = 4, DIR_ORI = 8,
+ DIR_ALL = 15, DIR_DNA = 15, DIR_PRO = 9
+ };
+
+ inline void add(const seq_data &s) {
+ // check to see if "name" already exists
+ // add to name if it does, otherwise push new sequence
+ std::pair<container::iterator,bool> res = cont.push_back(s);
+ using boost::lambda::_1;
+ if(!res.second)
+ cont.modify(res.first, _1 += s);
+ // remove gaps
+ cont.modify(res.first, _1 -= ss_gaps);
+ }
+
+ inline size_type size() const { return cont.size(); }
+
+ inline const seq_data& operator[](size_type sz) const {
+ return cont[sz];
+ }
+
+ inline void clear() {
+ cont.clear();
+ }
+
+ const container& db() const {
+ return cont;
+ }
+
+ template<typename InputIterator>
+ void rearrange(InputIterator it) {
+ cont.rearrange(it);
+ }
+
+ void transform(int dir) {
+ for(container::iterator it = cont.begin(); it != cont.end(); ++it)
+ cont.modify(it, bind(&seq_data::transform, _1, dir));
+ }
+
+ bool parse_file(const char *csfile, bool bappend=false, bool bi=false);
+ int unique_sort(int directions=DIR_ORI);
+
+ seq_db() : ss_gaps("-") { }
+ seq_db(const std::string &g) : ss_gaps(g) { }
+
+protected:
+ std::string ss_gaps;
+
+ container cont;
+};
+
+inline void seq_transform(std::string& dna, int dir) {
+ switch(dir&7) {
+ default:
+ case 0:
+ break;
+ case seq_db::DIR_COM:
+ std::transform(dna.begin(), dna.end(), dna.begin(), complement);
+ break;
+ case seq_db::DIR_RVC:
+ case (seq_db::DIR_REV|seq_db::DIR_COM):
+ std::transform(dna.begin(), dna.end(), dna.begin(), complement);
+ case seq_db::DIR_REV:
+ std::reverse(dna.begin(), dna.end());
+ break;
+ };
+}
+
+inline void seq_data::transform(int dir) {
+ seq_transform(dna, dir);
+}
+
+#endif
diff --git a/src/seqparser.h b/src/seqparser.h
new file mode 100644
index 0000000..86d4704
--- /dev/null
+++ b/src/seqparser.h
@@ -0,0 +1,150 @@
+/****************************************************************************
+ * Copyright (C) 2005-2010 Reed A. Cartwright, PhD <reed at scit.us> *
+ * *
+ * This program is free software: you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation, either version 3 of the License, or *
+ * (at your option) any later version. *
+ * *
+ * This program is distributed in the hope that it will be useful, *
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of *
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
+ * GNU General Public License for more details. *
+ * *
+ * You should have received a copy of the GNU General Public License *
+ * along with this program. If not, see <http://www.gnu.org/licenses/>. *
+ ****************************************************************************/
+
+#include <iostream>
+#include <stack>
+
+#include "seqdb.h"
+
+#include <boost/spirit/include/classic.hpp>
+#include <boost/spirit/include/classic_core.hpp>
+using namespace boost::spirit::classic;
+
+struct push_string
+{
+ push_string(std::stack<std::string>& st) : my_stack(st) { }
+
+ template<typename It> void operator() ( It first, It last) const
+ {
+ my_stack.push(std::string(first, last));
+ }
+
+ std::stack<std::string>& my_stack;
+};
+
+inline void trim_string(std::string &ss, const char *ch = " \n\t\r\v\f")
+{
+ std::string::size_type sz = ss.find_last_not_of(ch);
+ if(sz != std::string::npos) ss.erase(sz+1);
+ sz = ss.find_first_not_of(ch);
+ if(sz != std::string::npos) ss.erase(0, sz);
+}
+
+inline void sanitize_string(std::string &ss, const char *ch = " \n\t\r\v\f")
+{
+ std::string::size_type sz = ss.find_first_of(ch);
+ while(sz != std::string::npos)
+ {
+ ss.erase(sz, 1);
+ sz = ss.find_first_of(ch, sz);
+ }
+}
+
+struct pop_sequence
+{
+ pop_sequence(std::stack<std::string>& st, seq_db& db) : my_stack(st), rdb(db) { }
+
+ template<typename It> void operator() ( It first, It last) const
+ {
+ std::string seq(my_stack.top());
+ my_stack.pop();
+ std::string name(my_stack.top());
+ my_stack.pop();
+ trim_string(name);
+ sanitize_string(seq);
+
+ rdb.add(seq_data(name, seq));
+ }
+
+ std::stack<std::string>& my_stack;
+ seq_db &rdb;
+};
+
+struct add_sequence
+{
+ add_sequence(seq_db::size_type &pos, seq_db& db) : rpos(pos), rdb(db) { }
+
+ template<typename It> void operator() ( It first, It last) const
+ {
+ std::string seq(first,last);
+ sanitize_string(seq);
+ seq_db::size_type pos = rpos++ % rdb.size();
+ rdb.add(seq_data(rdb[pos].name, seq));
+ }
+
+ seq_db::size_type &rpos;
+ seq_db &rdb;
+};
+
+
+struct seq_grammar : public grammar<seq_grammar> {
+ seq_grammar(std::stack<std::string>& st, seq_db& db, seq_db::size_type &p) : string_stack(st), rdb(db), pos(p) {}
+
+ template <typename ScannerT> struct definition {
+ definition(seq_grammar const& self) : aln_sp("*:. ") {
+ self.pos = 0;
+
+ file_format = (*space_p) >> (phylip_format|aln_format|fasta_format)
+ >> (*space_p);
+ blank_line = *blank_p >> eol_p;
+
+ fasta_format = *fasta_seq;
+ fasta_seq = (fasta_seq_head >> fasta_seq_body)
+ [pop_sequence(self.string_stack, self.rdb)];
+ fasta_seq_head = ch_p('>')
+ >>(+(graph_p|blank_p))[push_string(self.string_stack)]
+ >> eol_p
+ ;
+ fasta_seq_body = (+(~ch_p('>')))[push_string(self.string_stack)];
+
+ aln_format = aln_head >> +aln_line;
+ aln_head = str_p("CLUSTAL") >> *(graph_p|blank_p) >> eol_p;
+ aln_line = (aln_seq | aln_special | *blank_p) >> eol_p;
+ aln_seq = (aln_seq_name >> +blank_p >> aln_seq_body >> *blank_p
+ >> *digit_p)
+ [pop_sequence(self.string_stack, self.rdb)];
+ aln_seq_name = (+graph_p)[push_string(self.string_stack)];
+ aln_seq_body = (+graph_p)[push_string(self.string_stack)];
+ aln_special = +aln_sp;
+
+ phylip_format = phylip_head >> phylip_block1
+ >> *(blank_line >> phylip_block);
+ phylip_head = uint_p >> +(blank_p) >> uint_p >> (*blank_p >> eol_p);
+ phylip_block1 = *(phylip_seq1 >> eol_p);
+ phylip_seq1 = (phylip_seq_name >> phylip_seq_seq)
+ [pop_sequence(self.string_stack, self.rdb)];
+ phylip_seq_name = (repeat_p(10)[graph_p|blank_p])[push_string(self.string_stack)];
+ phylip_seq_seq = (+(graph_p|blank_p))[push_string(self.string_stack)];
+ phylip_block = *(*blank_p >> phylip_seq >> eol_p);
+ phylip_seq = (+(graph_p|blank_p))[add_sequence(self.pos, self.rdb)];
+
+ }
+ chset<> aln_sp;
+ rule<ScannerT> file_format, blank_line;
+ rule<ScannerT> fasta_format, fasta_seq, fasta_seq_head, fasta_seq_body;
+ rule<ScannerT> aln_format, aln_head, aln_line, aln_seq, aln_seq_name,
+ aln_seq_body, aln_special;
+ rule<ScannerT> phylip_format, phylip_head, phylip_block1, phylip_block,
+ phylip_seq1, phylip_seq, phylip_seq_name, phylip_seq_seq;
+
+ rule<ScannerT> const& start() const { return file_format; }
+ };
+
+ std::stack<std::string> &string_stack;
+ seq_db::size_type &pos;
+ seq_db &rdb;
+};
diff --git a/src/sort.h b/src/sort.h
new file mode 100644
index 0000000..36e3ecd
--- /dev/null
+++ b/src/sort.h
@@ -0,0 +1,152 @@
+/****************************************************************************
+ * Copyright (C) 2009-2010 Reed A. Cartwright, PhD <reed at scit.us> *
+ * *
+ * This program is free software: you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation, either version 3 of the License, or *
+ * (at your option) any later version. *
+ * *
+ * This program is distributed in the hope that it will be useful, *
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of *
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
+ * GNU General Public License for more details. *
+ * *
+ * You should have received a copy of the GNU General Public License *
+ * along with this program. If not, see <http://www.gnu.org/licenses/>. *
+ ****************************************************************************/
+
+#pragma once
+#ifndef NGILA_SORT_H
+#define NGILA_SORT_H
+
+#include <algorithm>
+#include <vector>
+#include <string>
+#include <boost/functional/hash.hpp>
+
+enum {
+ DIR_ORI = 0,
+ DIR_REV = 1,
+ DIR_COM = 2,
+ DIR_RVC = 4,
+ DIR_ALL = 7
+};
+
+inline char dna_complement(char ch) {
+ static char dnac[] = "\0\x01\x02\x03\x04\x05\x06\a\b\t\n\x0b\f\r\x0e\x0f"
+ "\x10\x11\x12\x13\x14\x15\x16\x17\x18\x19\x1a\x1b\x1c\x1d\x1e\x1f"
+ " !\"#$%&'()*+,-./0123456789:;<=>?@"
+ "TVGHEFCDIJMLKNOPQYSAABWXRZ" "[\\]^_`"
+ "tvghefcdijmlknopqysaabwxrz" "{|}~\x7f";
+ return dnac[ch];
+}
+
+template<class It>
+inline void seq_transform(It first, It last, int dir) {
+ switch(dir&7) {
+ case DIR_REV:
+ std::reverse(first,last);
+ break;
+ case DIR_RVC:
+ case (DIR_COM|DIR_REV):
+ std::reverse(first, last);
+ case DIR_COM:
+ std::transform(first, last, first, dna_complement);
+ case DIR_ORI:
+ default:
+ break;
+ }
+}
+
+template<class Base>
+struct seq_sort_element {
+ seq_sort_element(const Base& s) : seq(s), direction(0) { }
+ Base seq;
+ std::size_t order;
+ std::size_t orig_order;
+ std::size_t direction;
+};
+
+template<class _Seq>
+class seq_sort : public std::vector< seq_sort_element<_Seq> > {
+public:
+ typedef seq_sort_element<_Seq> element_type;
+ typedef std::vector< seq_sort_element<_Seq> > base_type;
+ typedef _Seq value_type;
+
+ seq_sort() { }
+
+ template<class It>
+ seq_sort(It first, It last) : base_type(first, last) {
+ typename base_type::size_type o = 0;
+ for(typename base_type::iterator it = this->begin(); it != this->end(); ++it) {
+ it->order = it->orig_order = o++;
+ }
+ }
+
+ struct hash_data {
+ std::size_t hash_value;
+ std::size_t orig_order;
+ std::size_t direction;
+ hash_data(std::size_t a, std::size_t b, std::size_t c) :
+ hash_value(a), orig_order(b), direction(c) { }
+
+ bool operator<(const hash_data &h) const {
+ return hash_value < h.hash_value;
+ }
+ };
+
+ // Dirs is a flag that specifies which directions to try
+ int sort_seqs(int dirs) {
+ // create hashing object
+ boost::hash<value_type> hasher;
+ // store is a buffer of hashes, which will be used to determine
+ // unique ordering and direction
+ std::vector<hash_data> store;
+ std::size_t o = 0;
+ value_type ss;
+
+ // fill store with hashed values
+ for(typename base_type::iterator it = this->begin(); it != this->end(); ++it,++o) {
+ // original
+ store.push_back(hash_data(hasher(it->seq), o, DIR_ORI));
+ dirs &= 7;
+ if(dirs == DIR_ORI)
+ continue;
+ ss = it->seq;
+ if(dirs & DIR_COM == DIR_ORI) { // no complement sequence requested
+ seq_transform(ss.begin(), ss.end(), DIR_REV);
+ if(dirs & DIR_REV)
+ store.push_back(hash_data(hasher(ss), o, DIR_REV));
+ if(dirs & DIR_RVC) {
+ seq_transform(ss.begin(), ss.end(), DIR_COM);
+ store.push_back(hash_data(hasher(ss), o, DIR_RVC));
+ }
+ } else { // complement sequence is requested
+ seq_transform(ss.begin(), ss.end(), DIR_COM);
+ store.push_back(hash_data(hasher(ss), o, DIR_COM));
+ if(dirs != DIR_COM) { // other sequences requested
+ seq_transform(ss.begin(), ss.end(), DIR_REV);
+ if(dirs & DIR_RVC)
+ store.push_back(hash_data(hasher(ss), o, DIR_RVC));
+ if(dirs & DIR_REV) {
+ seq_transform(ss.begin(), ss.end(), DIR_COM);
+ store.push_back(hash_data(hasher(ss), o, DIR_REV));
+ }
+ }
+
+ }
+ }
+ // sort the hashed values
+ std::sort(store.begin(), store.end());
+
+ //lets see the look
+ for(typename std::vector<hash_data>::iterator it = store.begin(); it != store.end(); ++it) {
+ std::cerr << it->hash_value << " "
+ << it->orig_order << " "
+ << it->direction << std::endl;
+ }
+ }
+};
+
+#endif // NGILA_SORT_H
diff --git a/src/xm.h b/src/xm.h
new file mode 100644
index 0000000..29d3900
--- /dev/null
+++ b/src/xm.h
@@ -0,0 +1,76 @@
+/******************************************************************************
+ * Copyright (c) 2007-2009 Reed A. Cartwright, PhD <reed at scit.us> *
+ * *
+ * Permission is hereby granted, free of charge, to any person obtaining a *
+ * copy of this software and associated documentation files (the "Software"), *
+ * to deal in the Software without restriction, including without limitation *
+ * the rights to use, copy, modify, merge, publish, distribute, sublicense, *
+ * and/or sell copies of the Software, and to permit persons to whom the *
+ * Software is furnished to do so, subject to the following conditions: *
+ * *
+ * The above copyright notice and this permission notice shall be included in *
+ * all copies or substantial portions of the Software. *
+ * *
+ * 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 OR COPYRIGHT HOLDERS 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. *
+ ******************************************************************************/
+
+#include <boost/preprocessor.hpp>
+
+#ifndef XMACROS_HELPERS
+#define XMACROS_HELPERS
+
+/******************************************************************************
+ * X-Helpers List *
+ ******************************************************************************/
+
+// The _JS macro cats a seq 'seq' with separator 'sep'
+
+#define _JS_OP(s, data, elem) BOOST_PP_CAT(data, elem)
+
+#define _JS(sep, seq) BOOST_PP_IF( BOOST_PP_EQUAL(BOOST_PP_SEQ_SIZE(seq),1), \
+ _JS_1, _JS_2)(sep,seq)
+
+#define _JS_1(sep, seq) BOOST_PP_SEQ_HEAD(seq)
+
+#define _JS_2(sep, seq) BOOST_PP_SEQ_CAT(( BOOST_PP_SEQ_HEAD(seq) ) \
+ BOOST_PP_SEQ_TRANSFORM(_JS_OP, sep, BOOST_PP_SEQ_TAIL(seq)) \
+)
+
+// The _SS macro is similiar to _JS except that it stringizes everything
+
+#define _SS_OP(r, data, elem) data BOOST_PP_STRINGIZE(elem)
+
+#define _SS(sep, seq) BOOST_PP_STRINGIZE(BOOST_PP_SEQ_HEAD(seq)) \
+ BOOST_PP_SEQ_FOR_EACH(_SS_OP, sep, BOOST_PP_SEQ_TAIL(seq))
+
+// Output result if seq is defined
+
+#define _IFD(seq,res) BOOST_PP_EXPR_IF(BOOST_PP_GREATER(BOOST_PP_SEQ_SIZE((_)seq),1), res)
+
+#define _V(lname) _JS(_, lname)
+#define _S(lname) _SS("-", lname)
+#define _DL(a,b) a,b
+
+#else
+
+/******************************************************************************
+ * Cleanup *
+ ******************************************************************************/
+
+#undef XMACROS_HELPERS
+#undef _JS_OP
+#undef _JS
+#undef _SS
+#undef _IFD
+#undef _V
+#undef _S
+#undef _DL
+
+#endif
+
diff --git a/tests/CMakeLists.txt b/tests/CMakeLists.txt
new file mode 100644
index 0000000..e69de29
diff --git a/tests/test1.aln b/tests/test1.aln
new file mode 100644
index 0000000..fccd0ce
--- /dev/null
+++ b/tests/test1.aln
@@ -0,0 +1,36 @@
+A A----CTTGTCCGTTTCTGCAGGGAGTTAATATGCGTTTGGCTATGTACA-AATGTGAAA 55
+B ATCAGCTGGTCCGTTTCGGCAGGGAGCTAATATGCGTATGGCTATGTACATAATGTGAGA 60
+ * ** ********* ******** ********** ************ ******* *
+
+A ACGGTAGCCCTTTGTTAGGGTTTTTGTACGTCGGACGAACATTTGATGAGGTCATGGATG 115
+B ACGGTAGCCCTTTGTTAGGGTTTTTGTACGTCGGGCGAACATCTGAGGAGGTCATGGATC 120
+ ********************************** ******* *** ************
+
+A AATCCTATTCCAGAAGGTACGAATTCGGCTCACGTTGACCGCGTCGGGCACCTGCAATCG 175
+B AATCCTATGCCAGAAGGGAAGAATTCGGCTCACGTTGACCGCGCCGGGCACCGGCAATCG 180
+ ******** ******** * *********************** ******** *******
+
+A AATTATAATGAATTAACTTCCCACATCAGAGTTATCAAGCTTTCGGCTGGGTCTGCAA-- 233
+B AGTTATAATGAATTAACTCCCCACATCAGAGTTATCAAGCTTCCGGCTGGG-CTGCAAGT 239
+ * **************** *********************** ******** ******
+
+A ------------------------------------------------------------ 233
+B ACCACCTCCGCGACTCCACCTCCTCAGGAATAGAGCCCCCGTAACCAATTAAGAGCGGTT 299
+
+
+A ------------------TTGCGTTCCTTATCTGAGACTTGCTTTTGTCATCGCGAAATT 275
+B TCAGTTTGCCAGAGTCCGTTGCGTTCTTTATCTGAGACTTGGTTTTATTGTCGC--AATT 357
+ ******** ************** **** * **** ****
+
+A CCGGTCCTGGTTCCCGTGGAGTTCACTAGACTCTGGGGTATACCGGCGCTGTGCTTCACC 335
+B CCGGTTTTGGTTCCCGTGGAGTTCACTAGACTCTGGGG-------GCAGTGTGCGTCACC 410
+ ***** ******************************* ** ***** *****
+
+A CCTCAGTGGCTCGTAGGGTGTTGGTCCCCGGGACACCAGCTTAACAGAGCCATATATCAT 395
+B CCTCGGTGGCTCGTAGGGTGTTGGTCCTGGTGACACCAGCTTAAAAGAGGCCTATATCTT 470
+ **** ********************** * ************* **** * ****** *
+
+A GACTTCGTCGTCAACCCATTCAGTGTGAAA 425
+B TACTTCGTCGTCAACCCATCAAGTGTGAAA 500
+ ****************** *********
+
diff --git a/tests/test1.fsa b/tests/test1.fsa
new file mode 100644
index 0000000..823296a
--- /dev/null
+++ b/tests/test1.fsa
@@ -0,0 +1,22 @@
+>A
+A----CTTGTCCGTTTCTGCAGGGAGTTAATATGCGTTTGGCTATGTACA-AATGTGAAA
+ACGGTAGCCCTTTGTTAGGGTTTTTGTACGTCGGACGAACATTTGATGAGGTCATGGATG
+AATCCTATTCCAGAAGGTACGAATTCGGCTCACGTTGACCGCGTCGGGCACCTGCAATCG
+AATTATAATGAATTAACTTCCCACATCAGAGTTATCAAGCTTTCGGCTGGGTCTGCAA--
+------------------------------------------------------------
+------------------TTGCGTTCCTTATCTGAGACTTGCTTTTGTCATCGCGAAATT
+CCGGTCCTGGTTCCCGTGGAGTTCACTAGACTCTGGGGTATACCGGCGCTGTGCTTCACC
+CCTCAGTGGCTCGTAGGGTGTTGGTCCCCGGGACACCAGCTTAACAGAGCCATATATCAT
+GACTTCGTCGTCAACCCATTCAGTGTGAAA
+
+>B
+ATCAGCTGGTCCGTTTCGGCAGGGAGCTAATATGCGTATGGCTATGTACATAATGTGAGA
+ACGGTAGCCCTTTGTTAGGGTTTTTGTACGTCGGGCGAACATCTGAGGAGGTCATGGATC
+AATCCTATGCCAGAAGGGAAGAATTCGGCTCACGTTGACCGCGCCGGGCACCGGCAATCG
+AGTTATAATGAATTAACTCCCCACATCAGAGTTATCAAGCTTCCGGCTGGG-CTGCAAGT
+ACCACCTCCGCGACTCCACCTCCTCAGGAATAGAGCCCCCGTAACCAATTAAGAGCGGTT
+TCAGTTTGCCAGAGTCCGTTGCGTTCTTTATCTGAGACTTGGTTTTATTGTCGC--AATT
+CCGGTTTTGGTTCCCGTGGAGTTCACTAGACTCTGGG-------GGCAGTGTGCGTCACC
+CCTCGGTGGCTCGTAGGGTGTTGGTCCTGGTGACACCAGCTTAAAAGAGGCCTATATCTT
+TACTTCGTCGTCAACCCATCAAGTGTGAAA
+
diff --git a/tests/test1.sh b/tests/test1.sh
new file mode 100755
index 0000000..ba880e0
--- /dev/null
+++ b/tests/test1.sh
@@ -0,0 +1,10 @@
+#! /bin/sh
+
+: ${NGILA=ngila}
+: ${DIFF=diff}
+
+${NGILA} -a 2 -b 0.25 -c 0.5 ${top_srcdir}/tests/test1.fsa \
+ | tail -n +2 | ${DIFF} -b -B - ${top_srcdir}/tests/test1.aln
+result=$?
+
+exit $result
diff --git a/tools/aamodel.R b/tools/aamodel.R
new file mode 100644
index 0000000..5933f4f
--- /dev/null
+++ b/tools/aamodel.R
@@ -0,0 +1,25 @@
+# order of nucleotides in the file
+aa <- unlist(strsplit("ARNDCQEGHILKMFPSTWYV",""))
+
+# load table
+v <- scan("lg_LG.PAML.txt")
+u <- v[1:190]
+f <- v[191:210]
+s <- matrix(0,20,20)
+s[upper.tri(s)] <- u
+s <- t(s)
+s[upper.tri(s)] <- u
+
+# create subst matrix
+m <- t(s*f)
+diag(m) <- -apply(m,1,sum)
+m <- t(m*sqrt(f))/sqrt(f)
+
+# decompose
+em <- eigen(m,symmetric=TRUE)
+#em$values[abs(em$values) < 2*.Machine$double.eps] <- 0
+
+v <- c(sqrt(f),em$values,t(em$vectors))
+v <- matrix(sprintf("% .17e", v),ncol=4,byrow=T)
+cat(paste(apply(v,1,function(x) paste(x,collapse=", ")), collapse=",\n"))
+cat("\n")
diff --git a/tools/lg_LG.PAML.txt b/tools/lg_LG.PAML.txt
new file mode 100644
index 0000000..3adbdf6
--- /dev/null
+++ b/tools/lg_LG.PAML.txt
@@ -0,0 +1,22 @@
+
+0.425093
+0.276818 0.751878
+0.395144 0.123954 5.076149
+2.489084 0.534551 0.528768 0.062556
+0.969894 2.807908 1.695752 0.523386 0.084808
+1.038545 0.363970 0.541712 5.243870 0.003499 4.128591
+2.066040 0.390192 1.437645 0.844926 0.569265 0.267959 0.348847
+0.358858 2.426601 4.509238 0.927114 0.640543 4.813505 0.423881 0.311484
+0.149830 0.126991 0.191503 0.010690 0.320627 0.072854 0.044265 0.008705 0.108882
+0.395337 0.301848 0.068427 0.015076 0.594007 0.582457 0.069673 0.044261 0.366317 4.145067
+0.536518 6.326067 2.145078 0.282959 0.013266 3.234294 1.807177 0.296636 0.697264 0.159069 0.137500
+1.124035 0.484133 0.371004 0.025548 0.893680 1.672569 0.173735 0.139538 0.442472 4.273607 6.312358 0.656604
+0.253701 0.052722 0.089525 0.017416 1.105251 0.035855 0.018811 0.089586 0.682139 1.112727 2.592692 0.023918 1.798853
+1.177651 0.332533 0.161787 0.394456 0.075382 0.624294 0.419409 0.196961 0.508851 0.078281 0.249060 0.390322 0.099849 0.094464
+4.727182 0.858151 4.008358 1.240275 2.784478 1.223828 0.611973 1.739990 0.990012 0.064105 0.182287 0.748683 0.346960 0.361819 1.338132
+2.139501 0.578987 2.000679 0.425860 1.143480 1.080136 0.604545 0.129836 0.584262 1.033739 0.302936 1.136863 2.020366 0.165001 0.571468 6.472279
+0.180717 0.593607 0.045376 0.029890 0.670128 0.236199 0.077852 0.268491 0.597054 0.111660 0.619632 0.049906 0.696175 2.457121 0.095131 0.248862 0.140825
+0.218959 0.314440 0.612025 0.135107 1.165532 0.257336 0.120037 0.054679 5.306834 0.232523 0.299648 0.131932 0.481306 7.803902 0.089613 0.400547 0.245841 3.151815
+2.547870 0.170887 0.083688 0.037967 1.959291 0.210332 0.245034 0.076701 0.119013 10.649107 1.702745 0.185202 1.898718 0.654683 0.296501 0.098369 2.188158 0.189510 0.249313
+
+0.079066 0.055941 0.041977 0.053052 0.012937 0.040767 0.071586 0.057337 0.022355 0.062157 0.099081 0.064600 0.022951 0.042302 0.044040 0.061197 0.053287 0.012066 0.034155 0.069147
\ No newline at end of file
diff --git a/tools/logzeta.R b/tools/logzeta.R
new file mode 100644
index 0000000..4238b0c
--- /dev/null
+++ b/tools/logzeta.R
@@ -0,0 +1,30 @@
+library(gsl)
+
+x <- seq(1,10.1,0.1)
+y <- log((x-1)*zeta(x))
+y[1] <- 0
+f <- splinefun(x,y)
+z <- get("z", envir = environment(f))
+
+fmt <- "%.17e"
+
+m <- matrix(sprintf(fmt,z$y),ncol=4,byrow=T)
+cat("\tstatic double y[92] = {\n\t\t")
+cat(paste(apply(m,1,function(x) paste(x,collapse=", ")), collapse=",\n\t\t"))
+cat("\n\t};\n")
+
+m <- matrix(sprintf(fmt,z$b),ncol=4,byrow=T)
+cat("\tstatic double b[92] = {\n\t\t")
+cat(paste(apply(m,1,function(x) paste(x,collapse=", ")), collapse=",\n\t\t"))
+cat("\n\t};\n")
+
+m <- matrix(sprintf(fmt,z$c),ncol=4,byrow=T)
+cat("\tstatic double c[92] = {\n\t\t")
+cat(paste(apply(m,1,function(x) paste(x,collapse=", ")), collapse=",\n\t\t"))
+cat("\n\t};\n")
+
+m <- matrix(sprintf(fmt,z$d),ncol=4,byrow=T)
+cat("\tstatic double d[92] = {\n\t\t")
+cat(paste(apply(m,1,function(x) paste(x,collapse=", ")), collapse=",\n\t\t"))
+cat("\n\t};\n")
+
--
Alioth's /usr/local/bin/git-commit-notice on /srv/git.debian.org/git/debian-med/ngila.git
More information about the debian-med-commit
mailing list