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

# ------------------------------------------------------------------------------
# 
# Author: Gabriele Girelli
# Email: gigi.ga90@gmail.com
# Version: 1.4.4
# Date: 20170711
# Project: oligo characterization
# Description:  calculate melting temperature of a provide NA duplex
# 
# References:
#  [1] Freier et al, PNAS(83), 1986;
#  [2] Sugimoto et al, Biochemistry(34), 1995.
#  [3] Allawi & Santalucia, Biochemistry(36), 1997;
#  [4] SantaLucia, PNAS(95), 1998;
#  [5] Owczarzy et al, Biochemistry(43), 2004;
#  [6] Owczarzy et al, Biochemistry(47), 2008;
#  [7] McConaughy et al, Biochemistry(8), 1969;
#  [8] Wright et al, Appl. env. microbiol.(80), 2014.
# 
# ------------------------------------------------------------------------------



# DEPENDENCIES =================================================================

import argparse
import os
import sys

import oligo_melting as OligoMelt

# PARAMETERS ===================================================================

# Add script description
parser = argparse.ArgumentParser(
    description = '''
Calculate melting temperature of a DNA duplex at provided [oligo],
[Na+], [Mg2+]. Either provide an oligo sequence or a file with one oligo
per line (and use -F option). References:
 [1] Freier et al, PNAS(83), 1986;
 [2] Sugimoto et al, Biochemistry(34), 1995.
 [3] Allawi & Santalucia, Biochemistry(36), 1997;
 [4] SantaLucia, PNAS(95), 1998;
 [5] Owczarzy et al, Biochemistry(43), 2004;
 [6] Owczarzy et al, Biochemistry(47), 2008;
 [7] McConaughy et al, Biochemistry(8), 1969;
 [8] Wright et al, Appl. env. microbiol.(80), 2014.
''')

# Add mandatory arguments
parser.add_argument('seq', type = str, nargs = 1, help = '''
    DNA duplex sequence (one strand only) or path to a FASTA file (use with -F).
    ''')

# Add arguments with default value
parser.add_argument('-t', '--type', type = str, nargs = 1,
    help = '''Duplex type. Possible values: DNA:DNA (based on ref.3, default),
    RNA:RNA (based on ref.1), DNA:RNA (based on ref.2., given DNA sequence)
    or RNA:DNA (based on ref.2, given RNA sequence). The first nucleic acid type
    indicates the provided sequence.''',
    choices = OligoMelt.Duplex.NN_LABELS,
    default = [OligoMelt.Duplex.DEFAULT_NN_LABEL])
parser.add_argument('-o', '--oconc', metavar = "oligo_conc",
    type = float, nargs = 1,
    help = '''Oligonucleotide concentration [M].
    Default: %.2E M''' % OligoMelt.Duplex.DEFAULT_OLIGO_CONC,
    default = [OligoMelt.Duplex.DEFAULT_OLIGO_CONC])
parser.add_argument('-n', '--naconc', metavar = "na_conc",
    type = float, nargs = 1,
    help = '''Na+ concentration [M].
    Default: %.2E M''' % OligoMelt.Duplex.DEFAULT_NA_CONC,
    default = [OligoMelt.Duplex.DEFAULT_NA_CONC])
parser.add_argument('-m', '--mgconc', metavar = "mg_conc",
    type = float, nargs = 1,
    help = '''Mg2+ concentration [M]. Note: Mg2+ correction overwrites Na+
    correction. Default: %.2E M''' % OligoMelt.Duplex.DEFAULT_MG_CONC,
    default = [OligoMelt.Duplex.DEFAULT_MG_CONC])
parser.add_argument('-f', '--faconc', type = float, nargs = 1,
    metavar = 'fa_conc', help = '''Formamide concentration in perc.v,v.
    Default: %.2f''' % OligoMelt.Duplex.DEFAULT_FA_CONC,
    default = [OligoMelt.Duplex.DEFAULT_FA_CONC])
parser.add_argument('--fa-mode', type = str, nargs = 1,
    metavar = 'fa_mode', help = '''Mode of formamide correction. "mcconaughy"
    for classical -0.72%%FA correction from ref. 7, "wright" for single reaction
    model correction from ref.8 (default).''',
    choices = OligoMelt.Duplex.FA_MODE_LABELS,
    default = [OligoMelt.Duplex.DEFAULT_FA_MODE])
parser.add_argument('--fa-mvalue', type = str, nargs = 1, metavar = 'm',
    help = '''Specify the formamide m-value to be used with the wright
    correction model. Use either a single value "x" or two values "xL+y" where
    L is the probe length. Default: %s''' % OligoMelt.Duplex.DEFAULT_FA_MVALUE,
    default = [OligoMelt.Duplex.DEFAULT_FA_MVALUE])
parser.add_argument('--t-curve', type = float, nargs = 2,
    metavar = ('range', 'step'), help = '''Temperature range and step for
    melting curve generation. Use --make-curve option to generate the curve.
    Default: %.1f degC range and %.2f degC step.
    ''' % (OligoMelt.Duplex.DEFAULT_T_CURVE_RANGE,
        OligoMelt.Duplex.DEFAULT_T_CURVE_STEP), default = [
        OligoMelt.Duplex.DEFAULT_T_CURVE_RANGE,
        OligoMelt.Duplex.DEFAULT_T_CURVE_STEP])
parser.add_argument('--out-curve', type = str, nargs = 1, metavar = "outname",
    help = '''Path to output table containing tabulated curves.''',
    default = [OligoMelt.Duplex.DEFAULT_OUT_CURVE])

# Add flags
parser.add_argument('-C', '--celsius',
    dest = 'celsius', action = 'store_const',
    const = True, default = False,
    help = 'Output temperature in Celsius degrees. Default: Kelvin')
parser.add_argument('-v', '--verbose',
    dest = 'verbose', action = 'store_const',
    const = True, default = False,
    help = 'Verbose output.')

# Parse arguments
args = parser.parse_args()

# Assign to in-script variables ------------------------------------------------

# Input oligo
seq = args.seq[0]

# Concentrations
oligo_conc = args.oconc[0]
na_conc = args.naconc[0]
mg_conc = args.mgconc[0]

# Thermodynamic table
tt_mode = args.type[0]

# Temperature units
celsius = args.celsius

# File as input
use_file = os.path.isfile(seq)

# Verbose mode
is_verbose = args.verbose

# Melting curve
curve_range = args.t_curve[0]
curve_step = args.t_curve[1]
curve_outpath = args.out_curve[0]

# Formamide
fa_conc = args.faconc[0]
fa_mode = args.fa_mode[0]
fa_mval_s = args.fa_mvalue[0]

# Additional checks ------------------------------------------------------------

# Check proper curve step/range pair
if curve_step > curve_range / 2:
    sys.exit("!!!ERROR! Curve step must be smaller than curve range.")
curve_range -= curve_range % curve_step

# FUNCTIONS ====================================================================

# RUN ==========================================================================

# Build argument dictionary
data = {
    'oligo_conc' : oligo_conc,
    'na_conc' : na_conc,
    'mg_conc' : mg_conc,
    'fa_conc' : fa_conc,
    'fa_mode' : fa_mode,
    'fa_mval_s' : fa_mval_s,
    'tt_mode' : tt_mode,
    'celsius' : celsius,
    'is_verbose' : is_verbose,
    'curve_step' : curve_step,
    'curve_range' : curve_range,
    'curve_outpath' : curve_outpath,
    'silent' : False
}

# CALCULATE --------------------------------------------------------------------

if not use_file:
    # Single sequence case
    data['name'] = 'seq'
    data['seq'] = seq
    OligoMelt.Duplex.calc_tm(**data)
else:
    if not is_verbose:
        print("oligo_name\tdG\tdH\tdS\tTm\tSeq")

    # Input file case
    curr_head = ""
    curr_seq = ""
    with open(seq) as fin:
        for row in fin:
            if ">" == row[0]:
                if not 0 == len(curr_seq) and not 0 == len(curr_head):
                    # Calculate before moving to the next item
                    if " " in curr_head: curr_head = curr_head.split(" ")[0]
                    data['name'] = curr_head
                    data['seq'] = curr_seq
                    OligoMelt.Duplex.calc_tm(**data)

                curr_head = row[1:].strip()
                curr_seq = ""
            else:
                curr_seq += row.strip()

    # Calculate for last item
    if " " in curr_head: curr_head = curr_head.split(" ")[0]
    data['name'] = curr_head
    data['seq'] = curr_seq
    OligoMelt.Duplex.calc_tm(**data)

# END ==========================================================================

################################################################################
