uf-rc

#!/bin/sh
#
#  uf-rc - Reverse a/o complement the sequences in an unfasta stream.
#  Copyright (C) 2016  Marco van Zwetselaar <io@zwets.it>
#
#  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/>.
#
#  This utility is part of http://io.zwets.it/unfasta

# Nucleic Acids
#  A  adenosine          C  cytidine             G  guanine
#  T  thymidine          N  A/G/C/T (any)        U  uridine (EXCLUDED)
#  K  G/T (keto)         S  G/C (strong)         Y  T/C (pyrimidine) 
#  M  A/C (amino)        W  A/T (weak)           R  G/A (purine)        
#  B  G/T/C              D  G/A/T                H  A/C/T      
#  V  G/C/A              -  gap of indeterminate length (EXCLUDED)

NUCL_ALPHABET="AaCcGgTtNnKkMmSsWwYyRrBbVvDdHh"
COMP_ALPHABET="TtGgCcAaNnMmKkSsWwRrYyVvBbHhDd"

# Function to exit this script with an error message on stderr
err_exit() {
    echo "$(basename "$0"): $*" >&2
    exit 1
}

# Function to show usage information and exit
usage_exit() {
    echo "
Usage: $(basename $0) [OPTIONS] [FILE ...]

  Write to standard output the reverse complement of every sequence in each FILE.
  If no FILE is present or FILE is '-', read from standard input.

  OPTIONS
   -r, --reverse-only     Reverse only, do not complement.
   -c, --complement-only  Complement only, do not reverse.
   -m, --mark[-header]    Document operation by attaching '(uf:...)' to headers.

  Case is preserved.  Will also work for the degenerate nucleotide letters.
  E.g. complement of Y (pyrimidine, T/C) is R (purine, A/T).
" >&2
    exit ${1:-1}
}

# Defaults

REVERSE=1
COMPLEMENT=1
MARK=0

# Parse options

while [ $# -ne 0 -a "$(expr "$1" : '\(.\)..*')" = "-" ]; do
    case $1 in
    -r|--reverse-only)    COMPLEMENT=0 ;;
    -c|--complement-only) REVERSE=0 ;;
    -m|--mark*)           MARK=1 ;;
    --help)               usage_exit 0 ;;
    *) usage_exit ;;
    esac
    shift || usage_exit
done

# Do the work

gawk -bO -v P="$(basename "$0")" -v C=$COMPLEMENT -v R=$REVERSE -v M=$MARK '
    NR%2==1 { print $0 (M ? " (uf:" (R?"reverse":"") (C?"complement":"") ")" : "") }
    NR%2==0 {
        if (R) for (i = length($0); i >= 1; --i) print_maybe_comp()
        else   for (i = 1; i <= length($0); ++i) print_maybe_comp()
        printf "\n"
    }
    function print_maybe_comp(     c,p) {    # Haha, globals rule, no need to pass i
        c = substr ($0,i,1)
        if (C) {
            p = index ("'$NUCL_ALPHABET'", c)
            if (!p) {
                print P ": invalid character, cannot complement: " c
                exit 1
            }
            c = substr ("'$COMP_ALPHABET'", p, 1)
        }
        printf "%c", c
    }' "$@"

# vim: sts=4:sw=4:et:si:ai