#!/usr/bin/env python3
"""
Chapter 1: Python for Next Gen Biologist
Script: ch01_dna_basics.py
Description: DNA sequence validation, GC content calculation, nucleotide frequencies,
             RNA transcription, and reverse complement generation.
Usage:
    python ch01_dna_basics.py --seq ATGCGATCGATCGATCGATAGCTAGCTA
    python ch01_dna_basics.py --file sequence.fasta
"""

import argparse
import sys

def validate_dna(sequence: str) -> bool:
    """Check if the sequence contains only valid IUPAC DNA bases (A, C, G, T, N)."""
    valid_bases = set("ACGTUN")
    return set(sequence.upper()).issubset(valid_bases)

def calculate_gc_content(sequence: str) -> float:
    """Calculate GC percentage of a DNA sequence."""
    seq = sequence.upper()
    g_count = seq.count('G')
    c_count = seq.count('C')
    total_len = len(seq)
    if total_len == 0:
        return 0.0
    return round(((g_count + c_count) / total_len) * 100, 2)

def transcribe(dna_sequence: str) -> str:
    """Transcribe DNA sequence into messenger RNA (T -> U)."""
    return dna_sequence.upper().replace('T', 'U')

def reverse_complement(sequence: str) -> str:
    """Generate the reverse complement of a DNA sequence (5' -> 3')."""
    complement_map = str.maketrans("ACGTNacgtn", "TGCANtgcan")
    return sequence.translate(complement_map)[::-1]

def nucleotide_counts(sequence: str) -> dict:
    """Return a dictionary of nucleotide counts."""
    seq = sequence.upper()
    return {
        'A': seq.count('A'),
        'C': seq.count('C'),
        'G': seq.count('G'),
        'T': seq.count('T'),
        'N': seq.count('N')
    }

def main():
    parser = argparse.ArgumentParser(description="Chapter 1: DNA sequence metrics and reverse complement.")
    parser.add_argument("-s", "--seq", help="Input DNA sequence directly as string")
    parser.add_argument("-f", "--file", help="Path to file containing DNA sequence")
    args = parser.parse_args()

    seq = None
    if args.seq:
        seq = args.seq.strip()
    elif args.file:
        with open(args.file, 'r') as f:
            lines = [line.strip() for line in f if not line.startswith('>')]
            seq = "".join(lines)
    else:
        # Default demo sequence: Human Beta-Globin (HBB) 5' CDS fragment
        seq = "ATGGTGCACCTGACTCCTGAGGAGAAGTCTGCCGTTACTGCCCTGTGGGGCAAGGTGAACGTGGATGAAGTTGGTGGTGAGGCCCTGGGCAGGCTGCTGG"
        print("[!] No input provided. Using demo Human HBB CDS fragment.\n")

    if not validate_dna(seq):
        print("[ERROR] Sequence contains non-standard nucleotides.", file=sys.stderr)
        sys.exit(1)

    print("=== DNA Sequence Analysis Report ===")
    print(f"Sequence Length : {len(seq)} bp")
    print(f"GC Content      : {calculate_gc_content(seq)}%")
    counts = nucleotide_counts(seq)
    print(f"Base Frequencies: A={counts['A']} | C={counts['C']} | G={counts['G']} | T={counts['T']} | N={counts['N']}")
    print(f"mRNA Transcript : {transcribe(seq)[:50]}...")
    print(f"Rev-Complement  : {reverse_complement(seq)[:50]}...")

if __name__ == "__main__":
    main()