#!/usr/bin/env python3

# -*- coding: utf-8 -*-
#   Copyright (C) 2016 Mateusz Krzysztof Łącki and Michał Startek.
#
#   This file is part of MassTodon.
#
# License: see LICENCE file.

import argparse
import json
from os import listdir as ls, makedirs as mkdir
from os.path import join as pjoin

from masstodon.data.constants import eps, infinity
from masstodon.parse.blocked_fragments import parse as parse_blocked_fragments
from masstodon.parse.mods import parse_mods
from masstodon.masstodon import masstodon_single
from masstodon.parse.spectrum_file import parse as parse_spectrum_file

debug = False

p = argparse.ArgumentParser(description="Analyze one precursor with masstodon.")

p.add_argument( "-s", "--spectrum",
help="Path to the spectrum file: parseable extensions include '.txt' for ASCII, '.mzxml', or '.mzml' (case insensitive).",
                required=True)

p.add_argument( "-orbi", "--orbitrap",
                action='store_const',
                const=True,
                default=False,
                help="This is an Orbitrap spectrum.")

p.add_argument( "-f", "--fasta",
                type=str,
                help="Fasta sequence.",
                required=True)

p.add_argument( "-m", "--mods",
                nargs="+",
help = "Modifications of particular amino acids, like 11=C10H-2, 10C_carbo=Ag3, 1N=K2H-2")

p.add_argument( "-q", "--charge",
                type=int,
help="The initial charge of the precursor filtered out in MS1.",
                required=True,
                dest='q')

p.add_argument("-dq", "--distance_charges",
help="How many consecutive amino-acids are there per charge. Defaults to 5 (including the amino acid the charge is linked to).",
                type=int,
                default=5)

p.add_argument( "-tol", "--tolerance",
help="The tolerance in the m/z axis (relative or absolute). Valid expressions include 0.05Da, 0,06Th, 10ppm, 12.2ppm, 20.0mmu, 10mmu. ppm = parts per million, mmu = milli-mass unit, Da = daltons, Th = Thompsons (for Thor almighty's sake).",
                required=True)

p.add_argument( "-std", "--standard_deviations_count",
help="The number of standard deviations defining the tolerance in the coarse filtering stage.",
                type=float,
                default=3.0)

p.add_argument( "-fr", "--fragments",
                help="Only 'cz' accepted for now. Planning other fragmentation schemes, including inner fragments.",
                default='cz')

p.add_argument( "-bfr", "--blocked_fragments",
help="Fragments not included in the analysis, e.g. 'z5', or 'c13z21c12c53z1': no spaces between diffferent names.",
                type=parse_blocked_fragments,
                default='c0')

p.add_argument( "-ep", "--min_entry_prob",
help="The minimal probability an envelope has to scoop to be included in the deconvolution graph.",
                type=float,
                default=.8)

p.add_argument( "-ip", "--iso_prob",
help="The probability p in the p-optimal isotopic distribution. Defaults to 0.999.",
                type=float,
                default=.999)

p.add_argument( "-mi", "--min_intensity",
help="Experimental peaks with lower height will be trimmed. Defaults to a tincy-wincy nonzero number.",
                type=float,
                default=eps)

p.add_argument( "-dm", "--deconvolution_method",
help="nonnegative least squares (nnls) or quantile regression (quantile).",
                default="nnls",
                choices=["nnls", "quantile"])

p.add_argument( "--include_zero_intensities",
                action='store_const',
                const=True,
                default=False,
                help="Include zero intensities in the fitting of fragments. Don't do it for the Orbitrap, as it presumably cuts away some of the signal below an arbitrary detection limit.")

p.add_argument( "-t", "--target_folder",
help="Folder to write the results into. If not existing, will be created. Defaults to current folder.",
                default='.')

a = p.parse_args()

if debug:
    print(a.__dict__)

mz, i = parse_spectrum_file(a.spectrum)

m, timings = masstodon_single(
    mz, i, a.fasta, a.q,
    min_mz_diff=infinity,
    modifications=parse_mods(a.mods),
    orbitrap=a.orbitrap,
    threshold=a.tolerance,
    isotopic_coverage=a.iso_prob,
    min_prob=a.min_entry_prob,
    std_cnt=a.standard_deviations_count,
    get_timings=True,
    include_zero_intensities=a.include_zero_intensities)

t = a.target_folder

if t is not '.':
    mkdir(t, exist_ok=True)

print("Saving spectrum to where you started your python from.")

m.plotly(pjoin(t,"spectrum.html"),
         shape='rectangles',
         show=True)

m.write(t)
with open(pjoin(t, 'timings.json'), 'w') as h:
    json.dump(timings, h, indent=4)