<commit_before>#!/usr/bin/env python

# standard library imports
import os
import sys
import logging
import re
import hashlib

# 3rd party imports
import simplejson

# KBase imports
import biokbase.Transform.script_utils as script_utils


# transformation method that can be called if this module is imported
# Note the logger has different levels it could be run.  
# See: https://docs.python.org/2/library/logging.html#logging-levels
#
# The default level is set to INFO which includes everything except DEBUG
def transform(shock_service_url=None, handle_service_url=None, 
              output_file_name=None, input_directory=None, 
              working_directory=None, shock_id=None, handle_id=None, 
              input_mapping=None, fasta_reference_only=False, 
              level=logging.INFO, logger=None):
    """
    Converts FASTA file to KBaseGenomes.ContigSet json string.  
    Note the MD5 for the contig is generated by uppercasing the sequence.
    The ContigSet MD5 is generated by taking the MD5 of joining the sorted 
    list of individual contig's MD5s with a comma separator.

    Args:
        shock_service_url: A url for the KBase SHOCK service.
        handle_service_url: A url for the KBase Handle Service.
        output_file_name: A file name where the output JSON string should be stored.  
                          If the output file name is not specified the name will default 
                          to the name of the input file appended with '_contig_set'
        input_directory: The directory where files will be read from.
        working_directory: The directory the resulting json file will be written to.
        shock_id: Shock id for the fasta file if it already exists in shock
        handle_id: Handle id for the fasta file if it already exists as a handle
        input_mapping: JSON string mapping of input files to expected types.  
                       If you don't get this you need to scan the input 
                       directory and look for your files.
        fasta_reference_only: Creates a reference to the fasta file in Shock, but does not store the sequences in the workspace object.  Not recommended unless the fasta file is larger than 1GB. This is the default behavior for files that large.
        level: Logging level, defaults to logging.INFO.
        
    Returns:
        JSON file on disk that can be saved as a KBase workspace object.

    Authors:
        <NAME>, <NAME>
    """

    if logger is None:
        logger = script_utils.stderrlogger(__file__)
    
    logger.info("Starting conversion of FASTA to KBaseGenomes.ContigSet")
    token = os.environ.get('KB_AUTH_TOKEN')
        
    if input_mapping is None:
        logger.info("Scanning for FASTA files.")
    
        valid_extensions = [".fa",".fasta",".fna"]
    
        files = os.listdir(input_directory)
        fasta_files = [x for x in files if os.path.splitext(x)[-1] in valid_extensions]
            
        assert len(fasta_files) != 0
    
        logger.info("Found {0}".format(str(fasta_files)))

        input_file_name = os.path.join(input_directory,files[0])
    
        if len(fasta_files) > 1:
            logger.warning("Not sure how to handle multiple FASTA files in this context. Using {0}".format(input_file_name))
    else:
        input_file_name = os.path.join(os.path.join(input_directory, "FASTA.DNA.Assembly"), simplejson.loads(input_mapping)["FASTA.DNA.Assembly"])
        
                
    logger.info("Building Object.")
 
    if not os.path.isfile(input_file_name):
        raise Exception("The input file name {0} is not a file!".format(input_file_name))        

    if not os.path.isdir(args.working_directory):
        raise Exception("The working directory {0} is not a valid directory!".format(working_directory))        

    logger.debug(fasta_reference_only)

    # default if not too large
    contig_set_has_sequences = True 
    if fasta_reference_only:
        contig_set_has_sequences = False 

    fasta_filesize = os.stat(input_file_name).st_size
    if fasta_filesize > 1000000000:
        # Fasta file too large to save sequences into the ContigSet object.
        contigset_warn = """The FASTA input file seems to be too large. A ContigSet
                            object will be created without sequences, but will
                            contain a reference to the file."""
        logger.warning(contigset_warn) 
        contig_set_has_sequences = False 

    input_file_handle = open(input_file_name, 'r')
    
    fasta_header = None
    sequence_list = []
    fasta_dict = dict()
    first_header_found = False
    contig_set_md5_list = []
    # Pattern for replacing white space
    pattern = re.compile(r'\s+')
    sequence_exists = False
    
    valid_chars = "-AaCcGgTtUuWwSsMmKkRrYyBbDdHhVvNn"
    amino_acid_specific_characters = "PpLlIiFfQqEe" 

    for current_line in input_file_handle:
        if (current_line[0] == ">"):
            # found a header line
            # Wrap up previous fasta sequence
            if (not sequence_exists) and first_header_found:
                logger.error("There is no sequence related to FASTA record : {0}".format(fasta_header))        
                raise Exception("There is no sequence related to FASTA record : {0}".format(fasta_header))
            if not first_header_found:
                first_header_found = True
            else:
                # build up sequence and remove all white space
                total_sequence = ''.join(sequence_list)
                total_sequence = re.sub(pattern, '', total_sequence)
                if not total_sequence :
                    logger.error("There is no sequence related to FASTA record : {0}".format(fasta_header)) 
                    raise Exception("There is no sequence related to FASTA record : {0}".format(fasta_header))
                for character in total_sequence:
                    if character not in valid_chars:
                        if character in amino_acid_specific_characters:
                            raise Exception("This fasta file may have amino acids in it instead of the required nucleotides.")
                        raise Exception("This FASTA file has non nucleic acid characters : {0}".format(character))
                fasta_key = fasta_header.strip()
                contig_dict = dict() 
                contig_dict["id"] = fasta_key 
                contig_dict["length"] = len(total_sequence) 
                contig_dict["name"] = fasta_key 
                contig_dict["description"] = "Note MD5 is generated from uppercasing the sequence" 
                contig_md5 = hashlib.md5(total_sequence.upper()).hexdigest() 
                contig_dict["md5"] = contig_md5 
                contig_set_md5_list.append(contig_md5)
                 
                if contig_set_has_sequences: 
                    contig_dict["sequence"]= total_sequence
                else: 
                    contig_dict["sequence"]= ""
                
                fasta_dict[fasta_key] = contig_dict
               
                # get set up for next fasta sequence
                sequence_list = []
                sequence_exists = False
            
            fasta_header = current_line.replace('>','')
        else:
            sequence_list.append(current_line)
            sequence_exists = True

    input_file_handle.close()

    # wrap up last fasta sequence
    if (not sequence_exists) and first_header_found: 
        logger.error("There is no sequence related to FASTA record : {0}".format(fasta_header))        
        raise Exception("There is no sequence related to FASTA record : {0}".format(fasta_header)) 
    elif not first_header_found :
        logger.error("There are no contigs in this file") 
        raise Exception("There are no contigs in this file") 
    else: 
        # build up sequence and remove all white space      
        total_sequence = ''.join(sequence_list)
        total_sequence = re.sub(pattern, '', total_sequence)
        if not total_sequence :
            logger.error("There is no sequence related to FASTA record : {0}".format(fasta_header)) 
            raise Exception("There is no sequence related to FASTA record : {0}".format(fasta_header)) 

        for character in total_sequence: 
            if character not in valid_chars: 
                if character in amino_acid_specific_characters:
                    raise Exception("This fasta file may have amino acids in it instead of the required nucleotides.")
                raise Exception("This FASTA file has non nucleic acid characters : {0}".format(character))

        fasta_key = fasta_header.strip()
        contig_dict = dict()
        contig_dict["id"] = fasta_key 
        contig_dict["length"] = len(total_sequence)
        contig_dict["name"] = fasta_key
        contig_dict["description"] = "Note MD5 is generated from uppercasing the sequence" 
        contig_md5 = hashlib.md5(total_sequence.upper()).hexdigest()
        contig_dict["md5"]= contig_md5
        contig_set_md5_list.append(contig_md5)
        
        if contig_set_has_sequences: 
            contig_dict["sequence"] = total_sequence 
        else:
            contig_dict["sequence"]= ""
         
        fasta_dict[fasta_key] = contig_dict 


    if output_file_name is None:
        # default to input file name minus file extenstion adding "_contig_set" to the end
        base = os.path.basename(input_file_name)
        output_file_name = "{0}_contig_set.json".format(os.path.splitext(base)[0])
    
    contig_set_dict = dict()
    contig_set_dict["md5"] = hashlib.md5(",".join(sorted(contig_set_md5_list))).hexdigest()
    contig_set_dict["id"] = output_file_name
    contig_set_dict["name"] = output_file_name
    contig_set_dict["source"] = "KBase"
    contig_set_dict["source_id"] = os.path.basename(input_file_name) 
    contig_set_dict["contigs"] = [fasta_dict[x] for x in sorted(fasta_dict.keys())]

    if shock_id is None:
        shock_info = script_utils.upload_file_to_shock(logger, shock_service_url, input_file_name, token=token)
        shock_id = shock_info["id"]
    
    contig_set_dict["fasta_ref"] = shock_id

    # For future development if the type is updated to the handle_reference instead of a shock_reference

    # This generates the json for the object
    objectString = simplejson.dumps(contig_set_dict, sort_keys=True, indent=4)
    if sys.getsizeof(objectString) > 1000000000 :
        contig_set_dict["contigs"] = []
        objectString = simplejson.dumps(contig_set_dict, sort_keys=True, indent=4)
        logger.warning("The fasta file has a very large number of contigs thus resulting in an object being too large if " <commit_msg>Added check for file with no fasta sequences<commit_after>#!/usr/bin/env python

# standard library imports
import os
import sys
import logging
import re
import hashlib

# 3rd party imports
import simplejson

# KBase imports
import biokbase.Transform.script_utils as script_utils


# transformation method that can be called if this module is imported
# Note the logger has different levels it could be run.  
# See: https://docs.python.org/2/library/logging.html#logging-levels
#
# The default level is set to INFO which includes everything except DEBUG
def transform(shock_service_url=None, handle_service_url=None, 
              output_file_name=None, input_directory=None, 
              working_directory=None, shock_id=None, handle_id=None, 
              input_mapping=None, fasta_reference_only=False, 
              level=logging.INFO, logger=None):
    """
    Converts FASTA file to KBaseGenomes.ContigSet json string.  
    Note the MD5 for the contig is generated by uppercasing the sequence.
    The ContigSet MD5 is generated by taking the MD5 of joining the sorted 
    list of individual contig's MD5s with a comma separator.

    Args:
        shock_service_url: A url for the KBase SHOCK service.
        handle_service_url: A url for the KBase Handle Service.
        output_file_name: A file name where the output JSON string should be stored.  
                          If the output file name is not specified the name will default 
                          to the name of the input file appended with '_contig_set'
        input_directory: The directory where files will be read from.
        working_directory: The directory the resulting json file will be written to.
        shock_id: Shock id for the fasta file if it already exists in shock
        handle_id: Handle id for the fasta file if it already exists as a handle
        input_mapping: JSON string mapping of input files to expected types.  
                       If you don't get this you need to scan the input 
                       directory and look for your files.
        fasta_reference_only: Creates a reference to the fasta file in Shock, but does not store the sequences in the workspace object.  Not recommended unless the fasta file is larger than 1GB. This is the default behavior for files that large.
        level: Logging level, defaults to logging.INFO.
        
    Returns:
        JSON file on disk that can be saved as a KBase workspace object.

    Authors:
        <NAME>, <NAME>
    """

    if logger is None:
        logger = script_utils.stderrlogger(__file__)
    
    logger.info("Starting conversion of FASTA to KBaseGenomes.ContigSet")
    token = os.environ.get('KB_AUTH_TOKEN')
        
    if input_mapping is None:
        logger.info("Scanning for FASTA files.")
    
        valid_extensions = [".fa",".fasta",".fna"]
    
        files = os.listdir(input_directory)
        fasta_files = [x for x in files if os.path.splitext(x)[-1] in valid_extensions]
            
        assert len(fasta_files) != 0
    
        logger.info("Found {0}".format(str(fasta_files)))

        input_file_name = os.path.join(input_directory,files[0])
    
        if len(fasta_files) > 1:
            logger.warning("Not sure how to handle multiple FASTA files in this context. Using {0}".format(input_file_name))
    else:
        input_file_name = os.path.join(os.path.join(input_directory, "FASTA.DNA.Assembly"), simplejson.loads(input_mapping)["FASTA.DNA.Assembly"])
        
                
    logger.info("Building Object.")
 
    if not os.path.isfile(input_file_name):
        raise Exception("The input file name {0} is not a file!".format(input_file_name))        

    if not os.path.isdir(args.working_directory):
        raise Exception("The working directory {0} is not a valid directory!".format(working_directory))        

    logger.debug(fasta_reference_only)

    # default if not too large
    contig_set_has_sequences = True 
    if fasta_reference_only:
        contig_set_has_sequences = False 

    fasta_filesize = os.stat(input_file_name).st_size
    if fasta_filesize > 1000000000:
        # Fasta file too large to save sequences into the ContigSet object.
        contigset_warn = """The FASTA input file seems to be too large. A ContigSet
                            object will be created without sequences, but will
                            contain a reference to the file."""
        logger.warning(contigset_warn) 
        contig_set_has_sequences = False 

    input_file_handle = open(input_file_name, 'r')
    
    fasta_header = None
    sequence_list = []
    fasta_dict = dict()
    first_header_found = False
    contig_set_md5_list = []
    # Pattern for replacing white space
    pattern = re.compile(r'\s+')
    sequence_exists = False
    
    valid_chars = "-AaCcGgTtUuWwSsMmKkRrYyBbDdHhVvNn"
    amino_acid_specific_characters = "PpLlIiFfQqEe" 

    for current_line in input_file_handle:
        if (current_line[0] == ">"):
            # found a header line
            # Wrap up previous fasta sequence
            if (not sequence_exists) and first_header_found:
                logger.error("There is no sequence related to FASTA record : {0}".format(fasta_header))        
                raise Exception("There is no sequence related to FASTA record : {0}".format(fasta_header))
            if not first_header_found:
                first_header_found = True
            else:
                # build up sequence and remove all white space
                total_sequence = ''.join(sequence_list)
                total_sequence = re.sub(pattern, '', total_sequence)
                if not total_sequence :
                    logger.error("There is no sequence related to FASTA record : {0}".format(fasta_header)) 
                    raise Exception("There is no sequence related to FASTA record : {0}".format(fasta_header))
                for character in total_sequence:
                    if character not in valid_chars:
                        if character in amino_acid_specific_characters:
                            raise Exception("This fasta file may have amino acids in it instead of the required nucleotides.")
                        raise Exception("This FASTA file has non nucleic acid characters : {0}".format(character))
                fasta_key = fasta_header.strip()
                contig_dict = dict() 
                contig_dict["id"] = fasta_key 
                contig_dict["length"] = len(total_sequence) 
                contig_dict["name"] = fasta_key 
                contig_dict["description"] = "Note MD5 is generated from uppercasing the sequence" 
                contig_md5 = hashlib.md5(total_sequence.upper()).hexdigest() 
                contig_dict["md5"] = contig_md5 
                contig_set_md5_list.append(contig_md5)
                 
                if contig_set_has_sequences: 
                    contig_dict["sequence"]= total_sequence
                else: 
                    contig_dict["sequence"]= ""
                
                fasta_dict[fasta_key] = contig_dict
               
                # get set up for next fasta sequence
                sequence_list = []
                sequence_exists = False
            
            fasta_header = current_line.replace('>','')
        else:
            sequence_list.append(current_line)
            sequence_exists = True

    input_file_handle.close()

    # wrap up last fasta sequence
    if (not sequence_exists) and first_header_found: 
        logger.error("There is no sequence related to FASTA record : {0}".format(fasta_header))        
        raise Exception("There is no sequence related to FASTA record : {0}".format(fasta_header)) 
    elif not first_header_found :
        logger.error("There are no contigs in this file") 
        raise Exception("There are no contigs in this file") 
    else: 
        # build up sequence and remove all white space      
        total_sequence = ''.join(sequence_list)
        total_sequence = re.sub(pattern, '', total_sequence)
        if not total_sequence :
            logger.error("There is no sequence related to FASTA record : {0}".format(fasta_header)) 
            raise Exception("There is no sequence related to FASTA record : {0}".format(fasta_header)) 

        for character in total_sequence: 
            if character not in valid_chars: 
                if character in amino_acid_specific_characters:
                    raise Exception("This fasta file may have amino acids in it instead of the required nucleotides.")
                raise Exception("This FASTA file has non nucleic acid characters : {0}".format(character))

        fasta_key = fasta_header.strip()
        contig_dict = dict()
        contig_dict["id"] = fasta_key 
        contig_dict["length"] = len(total_sequence)
        contig_dict["name"] = fasta_key
        contig_dict["description"] = "Note MD5 is generated from uppercasing the sequence" 
        contig_md5 = hashlib.md5(total_sequence.upper()).hexdigest()
        contig_dict["md5"]= contig_md5
        contig_set_md5_list.append(contig_md5)
        
        if contig_set_has_sequences: 
            contig_dict["sequence"] = total_sequence 
        else:
            contig_dict["sequence"]= ""
         
        fasta_dict[fasta_key] = contig_dict 


    if output_file_name is None:
        # default to input file name minus file extenstion adding "_contig_set" to the end
        base = os.path.basename(input_file_name)
        output_file_name = "{0}_contig_set.json".format(os.path.splitext(base)[0])
    
    contig_set_dict = dict()
    contig_set_dict["md5"] = hashlib.md5(",".join(sorted(contig_set_md5_list))).hexdigest()
    contig_set_dict["id"] = output_file_name
    contig_set_dict["name"] = output_file_name
    contig_set_dict["source"] = "KBase"
    contig_set_dict["source_id"] = os.path.basename(input_file_name) 
    contig_set_dict["contigs"] = [fasta_dict[x] for x in sorted(fasta_dict.keys())]

    if shock_id is None:
        shock_info = script_utils.upload_file_to_shock(logger, shock_service_url, input_file_name, token=token)
        shock_id = shock_info["id"]
    
    contig_set_dict["fasta_ref"] = shock_id

    # For future development if the type is updated to the handle_reference instead of a shock_reference

    # This generates the json for the object
    objectString = simplejson.dumps(contig_set_dict, sort_keys=True, indent=4)
    if len(contig_set_dict["contigs"]) == 0:
        raise Exception("There appears to be no FASTA DNA Sequences in the input file.") 
    if sys.getsizeof(objectString) > 1000000000 :
        contig_set_dict["contigs"] = []
        objectString = simplejson.dumps(contig_set_dict, sort_keys=True, indent=4)
        logger.warning("The fasta file has a very large number of contigs thus resulting in an object being too large if " 
{'fraction_non_alphanumeric': 0.07465582555668483, 'fraction_numerical': 0.006329113924050633, 'mean_word_length': 2.456330831869711, 'pattern_counts': {'":': 0, '<': 7, '<?xml version=': 0, '>': 17, 'https://': 2, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 4, 'word_list_counts': {'cursed_substrings.json': 4, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'dir': 'git-commits-cleaned', 'id': '834251', 'n_tokens_mistral': 6291, 'n_tokens_neox': 5923, 'n_words': 1821}