def dqueries2queriessam(cfg,dqueries):    
    """
    Aligns queries to genome and gets SAM file
    step#1

    :param cfg: configuration dict
    :param dqueries: dataframe of queries
    """
    datatmpd=cfg['datatmpd']
    dqueries=set_index(dqueries,'query id')
    queryls=dqueries.loc[:,'query sequence'].apply(len).unique()
    for queryl in queryls:
        logging.debug(f"now aligning queries of length {queryl}")
        queriesfap = f'{datatmpd}/01_queries_queryl{queryl:02}.fa'
        logging.info(basename(queriesfap))
        if not exists(queriesfap) or cfg['force']:
            with open(queriesfap,'w') as f:
                for gi in dqueries.index:
                    f.write('>{}\n{}\n'.format(gi.replace(' ','_'),dqueries.loc[gi,'query sequence']))
        ## BWA alignment command is adapted from cripror 
        ## https://github.com/rraadd88/crisporWebsite/blob/master/crispor.py
        # BWA allow up to X mismatches
        # maximum number of occurences in the genome to get flagged as repeats. 
        # This is used in bwa samse, when converting the sam file
        # and for warnings in the table output.
        MAXOCC = 60000

        # the BWA queue size is 2M by default. We derive the queue size from MAXOCC
        MFAC = 2000000/MAXOCC

        genomep=cfg['genomep']
        genomed = dirname(genomep) # make var local, see below
        genomegffp=cfg['genomegffp']

        # increase MAXOCC if there is only a single query, but only in CGI mode
        bwaM = MFAC*MAXOCC # -m is queue size in bwa
        queriessap = f'{datatmpd}/01_queries_queryl{queryl:02}.sa'
        logging.info(basename(queriessap))
        if not exists(queriessap) or cfg['force']:
            cmd=f"{cfg['bwa']} aln -t 1 -o 0 -m {bwaM} -n {cfg['mismatches_max']} -k {cfg['mismatches_max']} -N -l {queryl} {genomep} {queriesfap} > {queriessap} 2> {queriessap}.log"
            runbashcmd(cmd)

        queriessamp = f'{datatmpd}/01_queries_queryl{queryl:02}.sam'
        logging.info(basename(queriessamp))        
        if not exists(queriessamp) or cfg['force']:
            cmd=f"{cfg['bwa']} samse -n {MAXOCC} {genomep} {queriessap} {queriesfap} > {queriessamp} 2> {queriessamp}.log"
            runbashcmd(cmd)
    return cfg
def dqueries2queriessam ( cfg , dqueries ) : datatmpd = cfg [ 'datatmpd' ] dqueries = set_index ( dqueries , 'query id' ) queryls = dqueries . loc [ : , 'query sequence' ] . apply ( len ) . unique ( ) for queryl in queryls : logging . debug ( f"now aligning queries of length {queryl}" ) queriesfap = f'{datatmpd}/01_queries_queryl{queryl:02}.fa' logging . info ( basename ( queriesfap ) ) if not exists ( queriesfap ) or cfg [ 'force' ] : with open ( queriesfap , 'w' ) as f : for gi in dqueries . index : f . write ( '>{}\n{}\n' . format ( gi . replace ( ' ' , '_' ) , dqueries . loc [ gi , 'query sequence' ] ) ) MAXOCC = 60000 MFAC = 2000000 / MAXOCC genomep = cfg [ 'genomep' ] genomed = dirname ( genomep ) genomegffp = cfg [ 'genomegffp' ] bwaM = MFAC * MAXOCC queriessap = f'{datatmpd}/01_queries_queryl{queryl:02}.sa' logging . info ( basename ( queriessap ) ) if not exists ( queriessap ) or cfg [ 'force' ] : cmd = f"{cfg['bwa']} aln -t 1 -o 0 -m {bwaM} -n {cfg['mismatches_max']} -k {cfg['mismatches_max']} -N -l {queryl} {genomep} {queriesfap} > {queriessap} 2> {queriessap}.log" runbashcmd ( cmd ) queriessamp = f'{datatmpd}/01_queries_queryl{queryl:02}.sam' logging . info ( basename ( queriessamp ) ) if not exists ( queriessamp ) or cfg [ 'force' ] : cmd = f"{cfg['bwa']} samse -n {MAXOCC} {genomep} {queriessap} {queriesfap} > {queriessamp} 2> {queriessamp}.log" runbashcmd ( cmd ) return cfg
Aligns queries to genome and gets SAM file
    step#1

    :param cfg: configuration dict
    :param dqueries: dataframe of queries