#define PROG_NAME "dnasample"
#define PROG_DESC "extract sample segments of a DNA sequence"
#define PROG_VERS "1.1"
/* Last edited on 2024-12-21 11:20:52 by stolfi */

/* Copyright © 2003 by the the Fluminense Federal University (UFF).
** See the copyright, authorship, and warranty notice at end of file.
** Last edited on 2006-04-07 18:35:18 by stolfi
*/

#define PROG_HELP \
  PROG_NAME " INICIO NUM_SEQS TAM_SEQ PASSO \\\n" \
  "  < ARQ_FASTA \\\n" \
  "  > ARQ_SAIDA"

/*
** Este programa lê da entrada padrão {stdin} uma biosseqüência (DNA
** ou proteína), formada por letras A-Z; extrai da mesma um certo
** número {NUM_SEQS} de trechos; e escreve esses trechos na saída
** padrão {stdout}.
**
** Cada trecho tem um número fixo {TAM_SEQ} de letras. O usuário
** especifica o deslocamento nominal {INICIO} do primeiro trecho
** relativa ao início da entrada, e o deslocamento nominal {PASSO}
** entre os inícios de dois trechos consecutivos.
** 
** Em particular, se {PASSO} for menor que {TAM_SEQ}, cada trecho
** começará com as últimas {TAM_SEQ - PASSO} letras do trecho
** anterior. Se {PASSO} for maior que {TAM_SEQ}, os trechos serão
** sempre disjuntos. Neste caso, o início de cada trecho será
** aleatoriamente deslocado de sua posição nominal, mantendo-se
** a média {PASSO}.
**
** O arquivo de entrada deve estar no formato FASTA. O programa ignora
** brancos, fins-de-linha, e linhas que começam com ">". O arquivo de
** saída contém apenas as letras dos trechos, formatadas 70 por linha,
** com uma linha em branco ao fim de cada trecho.
** 
*/

#include <string.h>
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <jsrandom.h>

#define MAX_BASES 10000000

/* IMPLEMENTAÇÕES */

int main (int argc, char **argv)
{
  if (argc != 5)
    {
      fprintf(stderr, "Parametros errados!\n");
      fprintf(stderr, "%s versão %s, uso:\n%s\n", PROG_NAME, PROG_VERS, PROG_HELP);
      return(1);
    }

  int inicio = atoi(argv[1]);
  int numTrechos = atoi(argv[2]);
  int tamTrecho = atoi(argv[3]);
  int passo = atoi(argv[4]);
  fprintf(stderr, 
    "Gerando %d trechos com %d bases cada e passo %d...\n",
    numTrechos, tamTrecho, passo);

  /* Lê a cadeia de entrada, coloca em {base[0..tamEntrada-1]}: */
  char base[MAX_BASES]; /* Cadeia de entrada, limpa. */
  int tamEntrada = 0;   /* Numero de bases lidas. */    
  int nLinha;           /* Número da linha corrente. */
  int nColuna;          /* Número do caracter dentro da linha corrente. */
  int ignoraLinha;      /* Se 1, a linha corrente deve ser ignorada. */
  int letra;
  
  letra = getc(stdin);
  nLinha = 1; nColuna = 1; ignoraLinha = 0;
  while (letra != EOF)
    { 
      if (letra == '\n')
        { /* Fim de linha, prepara para começar outra: */
          nLinha++; nColuna = 0; ignoraLinha = 0;
        }
      else if (ignoraLinha)
        { /* Cabeçalho/comentário, ignora: */
        }
      else if ((letra == '>') && (nColuna == 1))
        { /* Começo de cabeçalho/comentário, ignora resto da linha: */
          ignoraLinha = 1;
        }
      else if (letra == ' ')
        { /* Branco, ignora. */
        }
      else if ((letra >= 'A') && (letra <= 'Z'))
        { 
          if (tamEntrada >= MAX_BASES) 
            { fprintf(stderr, 
                "Tamanho maximo da cadeia (%d) excedido\n", 
                MAX_BASES);
              return(1);
            }
          base[tamEntrada] = letra; tamEntrada++;
        }
      else
        { 
          fprintf(stderr, 
            "Linha %d: Caracter inválido (%c)\n",
            nLinha, letra);
          return(1);
        }
      letra = getc(stdin);
      nColuna++;
    }
      
  /* Deslocamento aleatório mínimo e máximo de cada amostra. */
  int kMin, kMax; 
  if (passo <= tamTrecho) 
    { 
      /* O deslocamento entre os trechos amostrados é fixo: */
      kMin = kMax = 0;
    }
  else
    { 
      /* O deslocamento é semi-aleatório, evitando contato: */
      kMin = (passo - tamTrecho)/4;
      kMax = passo - kMin - tamTrecho;
    }
    
  /* Extrai os trechos: */
  int p, i;
  for (p = 0; p < numTrechos; p++)
    {
      /* Escolhe a posição de partida: */
      int k = int32_abrandom(kMin, kMax);
      int inicioTrecho = inicio + p * passo + k;
      int fimTrecho = inicioTrecho + tamTrecho - 1;
      fprintf(stderr, "trecho %d = base[%d .. %d]\n", 
        p, inicioTrecho, fimTrecho);
      if (fimTrecho >= tamEntrada)
        { 
          fprintf(stderr, 
            "Trecho extrapola a cadeia de entrada\n");
          return(1);
        }
      for (i = 0; i < tamTrecho; i++)
	{ 
          if ((i > 0) && (i % 70 == 0)) { fputc('\n', stdout); }
          fputc(base[inicioTrecho + i], stdout);
        }
      fputc('\n', stdout);
      fputc('\n', stdout);
    }
  fclose(stdout);
  fclose(stdin);
  return 0;
}

/* COPYRIGHT, AUTHORSHIP, AND WARRANTY NOTICE:
**
**   Copyright © 2003 by the the Fluminense Federal University (UFF).
**   
** Created 2003 by Luciana Pessôa, IC-UFF.
** Heavily modified 2004-2005 by Jorge Stolfi, IC-UNICAMP.
**
** Permission to use, copy, modify, and redistribute this software and
** its documentation for any purpose and without fee is hereby
** granted, provided that: (1) the copyright notice at the top of this
** file and this copyright, authorship, and warranty notice is retained
** in all derived source files and documentation; (2) no executable
** code derived from this file is published or distributed without the
** corresponding source code; and (3) these same rights are granted to
** any recipient of such code, under the same conditions.
** This software is provided "as is", WITHOUT ANY EXPLICIT OR IMPLICIT
** WARRANTIES, not even the implied warranties of merchantibility and
** fitness for a particular purpose. END OF NOTICE.
*/