Commit 6d8bb449 by Eric Coissac

delete old version

git-svn-id: https://www.grenoble.prabi.fr/svn/LECASofts/ecoPrimers/trunk@179 60f365c0-8329-0410-b2a4-ec073aeeaa1d
parent f9410e5a
EXEC=ecoPrimer
PRIMER_SRC= ecoprimer.c
PRIMER_OBJ= $(patsubst %.c,%.o,$(PRIMER_SRC))
SRCS= $(PRIMER_SRC)
LIB= -lecoPCR -lapat -lKMRK -lz -lm
LIBFILE= libapat/libapat.a \
libecoPCR/libecoPCR.a \
libKMRK/libKMRK.a
include global.mk
all: $(EXEC)
########
#
# ecoPrimer compilation
#
########
# executable compilation and link
ecoPrimer: $(PRIMER_OBJ) $(LIBFILE)
$(CC) $(LDFLAGS) -o $@ $< $(LIBPATH) $(LIB)
########
#
# library compilation
#
########
libapat/libapat.a:
$(MAKE) -C libapat
libecoPCR/libecoPCR.a:
$(MAKE) -C libecoPCR
libKMRK/libKMRK.a:
$(MAKE) -C libKMRK
########
#
# project management
#
########
clean:
rm -f *.o
rm -f $(EXEC)
$(MAKE) -C libapat clean
$(MAKE) -C libecoPCR clean
$(MAKE) -C libKMRK clean
\ No newline at end of file
MACHINE=MAC_OS_X
LIBPATH= -Llibapat -LlibecoPCR -LlibKMRK
MAKEDEPEND = gcc -D$(MACHINE) -M $(CPPFLAGS) -o $*.d $<
CC=gcc
CFLAGS= -W -Wall -O2 -g
default: all
%.o: %.c
$(CC) -D$(MACHINE) $(CFLAGS) -c -o $@ $<
%.P : %.c
$(MAKEDEPEND)
@sed 's/\($*\)\.o[ :]*/\1.o $@ : /g' < $*.d > $@; \
rm -f $*.d; [ -s $@ ] || rm -f $@
include $(SRCS:.c=.P)
\ No newline at end of file
#ifndef KMRK_Seeds_h
#define KMRK_Seeds_h
/********************************************
********************************************
**
** Declaration of struct
**
********************************************
********************************************/
#include <stdint.h>
#include <stdio.h>
#include "KMRK.h"
#define KMRK_SORT_SEEDS_FUNC_PTR(name) int32_t (*name)(SmallSeed_type*, \
SmallSeed_type*)
#define KMRK_DELTA_SEEDS_FUNC_PTR(name) int32_t (*name)(SmallSeed_type*)
/********************************************
********************************************
**
** Declaration of public functions
**
********************************************
********************************************/
AllSeeds_type *KMRK_allocSeeds(AllSeeds_type *AllSeeds,
int32_t size,
int8_t opt_dir,
int8_t opt_inv);
void KMRK_SetMultipleLenDirSeeds(SmallSeed_type* seeds,
int32_t nseeds,
int32_t wordSize,
AllSeeds_type *PtrAllSeeds);
void KMRK_freeSeeds(AllSeeds_type *AllSeeds);
void KMRK_compactSeeds(AllSeeds_type *AllSeeds);
void KMRK_pushSeed(AllSeeds_type *AllSeeds,
int32_t pos1,
int32_t pos2,
int32_t length,
int8_t dir);
AllSeeds_type* KMRK_enumerateDirectCouple(AllSeeds_type* Seeds,
int32_t expected,
int32_t wordSize,
vn_type* stack,
int32_t seq);
AllSeeds_type* KMRK_enumerateInvertedCouple(AllSeeds_type* Seeds,
int32_t expected,
int32_t wordSize,
vn_type* stack);
AllSeeds_type* KMRK_enumerateInterCouple(AllSeeds_type* Seeds,
int32_t seq1,
int32_t seq2,
int32_t expected,
int32_t wordSize,
vn_type* stack);
AllSeeds_type* KMRK_enumerateInterInvertedCouple(AllSeeds_type* Seeds,
int32_t seq2,
int32_t expected,
int32_t wordSize,
vn_type* stack);
/**
* Compare two seeds and return an integer less than, equal to or greater
* than zero considering the relative order of the two seeds. This
* version take into account only pos1 and pos2 of seeds without taking
* account of the sequences or of the relative direction
*
* @param s1 pointer to seed one
* @param s2 pointer to seed two
*
* @return a integer less than, equal to or greater than zero
*/
int32_t KMRK_cmpSeedsPos(SmallSeed_type *s1, SmallSeed_type *s2);
int32_t KMRK_cmpDeltaSeedsPos(SmallSeed_type *s1, SmallSeed_type *s2);
int32_t KMRK_cmpDeltaInvSeedsPos(SmallSeed_type *s1, SmallSeed_type *s2);
void KMRK_sortSeeds(SmallSeed_type* seeds,
int32_t nseeds,
KMRK_SORT_SEEDS_FUNC_PTR(compare));
AllSeeds_type* KMRK_get_seeds(char **seq,
int32_t SimpleSeqLen,
int16_t Lmin,
int8_t opt_dir,
int8_t opt_inv,
int8_t opt_verbose,
masked_area_table_t *mask);
AllSeeds_type* KMRK_get_seeds_2seqs(char **seq1,
char **seq2,
int32_t size1,
int32_t size2,
int16_t Lmin,
int8_t opt_dir,
int8_t opt_inv,
int8_t opt_verbose,
masked_area_table_t *mask);
/**
* Order an array of seeds by pos1,pos2
*
* @param seeds pointer to an array of Seed_type object to sort
* @param count count of element in the array
*/
void KMRK_sortSeedsByPos(Seed_type* seeds, int32_t count);
#endif /* KMRK_Seeds_h */
/*
* KMRK_mask.c
* repseek
*
* Created by Eric Coissac on 04/12/04.
* Copyright 2004 __MyCompanyName__. All rights reserved.
*
*/
#include "KMRK_mask.h"
#include <stdio.h>
#include <stdlib.h>
#include "memory.h"
#define MASKED_AREA_TABLE_SIZE(seqcount) (sizeof(masked_area_table_t) + (sizeof(masked_area_list_t*) * ((seqcount)-1)))
#define MASKED_AREA_LIST_SIZE(areacount) (sizeof(masked_area_list_t) + (sizeof(masked_area_t) * ((areacount)-1)))
#define AREA_COUNT_INIT (1000)
static masked_area_table_t *new_masked_area_table(int32_t seqcount, int32_t areacount);
static masked_area_list_t *new_masked_area_list(int32_t areacount);
static masked_area_list_t *realloc_masked_area_list(masked_area_list_t *list,int32_t areacount);
static int32_t push_area(masked_area_table_t* table,int32_t sequence,int32_t begin,int32_t end);
static void sort_area_table(masked_area_table_t* table);
static int32_t compare_area(const masked_area_t* v1,const masked_area_t* v2);
static int32_t search_area(const masked_area_t* v1,const masked_area_t* v2);
static masked_area_list_t *strip_area_list(masked_area_list_t* list);
static void strip_area_table(masked_area_table_t* table);
static masked_area_list_t *new_masked_area_list(int32_t areacount)
{
masked_area_list_t *list;
list = MyCalloc(1,MASKED_AREA_LIST_SIZE(areacount),"Not enougth memory for mask table");
list->reserved=areacount;
return list;
}
static masked_area_list_t *realloc_masked_area_list(masked_area_list_t *list,int32_t areacount)
{
list = MyRealloc(list,
MASKED_AREA_LIST_SIZE(areacount),
MASKED_AREA_LIST_SIZE(list->reserved),
"Not enougth memory for mask table");
list->reserved=areacount;
return list;
}
static masked_area_table_t *new_masked_area_table(int32_t seqcount, int32_t areacount)
{
masked_area_table_t *table;
int32_t i;
table = MyCalloc(1,MASKED_AREA_TABLE_SIZE(seqcount),"Not enougth memory for mask table");
table->seqcount=seqcount;
for (i=0;i<seqcount;i++)
table->sequence[i]=new_masked_area_list(areacount);
return table;
}
static int32_t push_area(masked_area_table_t* table,int32_t sequence,int32_t begin,int32_t end)
{
masked_area_list_t * list;
if (sequence >= table->seqcount)
return -1;
list = table->sequence[sequence];
if (list->reserved == list->count)
{
list = realloc_masked_area_list(list,list->reserved*2);
table->sequence[sequence]=list;
}
list->area[list->count].begin=begin;
list->area[list->count].end=end;
list->count++;
table->total++;
return table->total;
}
static int32_t compare_area(const masked_area_t* v1,const masked_area_t* v2)
{
return v1->begin - v2->begin;
}
static void sort_area_table(masked_area_table_t* table)
{
int32_t i;
for (i=0; i<table->seqcount;i++)
{
qsort(table->sequence[i]->area,
table->sequence[i]->count,
sizeof(masked_area_t),
(int (*)(const void *, const void *))compare_area);
}
}
static masked_area_list_t *strip_area_list(masked_area_list_t* list)
{
int32_t i;
int32_t j;
int32_t count;
int32_t newcount;
count = list->count;
newcount=count;
for (i=1;i<count;i++)
{
/* fprintf(stderr,"\n%d->%d %d->%d ==>",list->area[i-1].begin,list->area[i-1].end,list->area[i].begin,list->area[i].end); */
if ((list->area[i].begin-1) <= list->area[i-1].end)
{
/* fprintf(stderr," joined"); */
list->area[i].begin=list->area[i-1].begin;
list->area[i-1].begin=-1;
newcount--;
}
}
list->count=newcount;
for (i=0,j=0;i<count;i++)
{
if (list->area[i].begin>=0)
{
if (i!=j)
list->area[j]=list->area[i];
j++;
}
}
return realloc_masked_area_list(list,newcount);
}
static void strip_area_table(masked_area_table_t* table)
{
int32_t seq;
int32_t oldcount;
masked_area_list_t* list;
sort_area_table(table);
for (seq=0; seq < table->seqcount;seq++)
{
list = table->sequence[seq];
oldcount = list->count;
table->sequence[seq]=strip_area_list(list);
table->total-=oldcount - table->sequence[seq]->count;
}
}
static int32_t search_area(const masked_area_t* v1,const masked_area_t* v2)
{
int32_t pos;
pos = v1->begin;
if (pos < v2->begin)
return -1;
if (pos > v2->end)
return 1;
return 0;
}
masked_area_table_t *KMRK_ReadMaskArea(char *areafile,int32_t seqcount,int32_t complement)
{
FILE* area;
char buffer[1000];
char *ok;
int32_t begin;
int32_t end;
int32_t sequence;
int32_t column;
int32_t linecount;
masked_area_table_t *table;
if (complement > 0)
seqcount++;
else
complement=0;
area = fopen(areafile,"r");
linecount=0;
table=new_masked_area_table(seqcount,AREA_COUNT_INIT);
do {
linecount++;
ok = fgets(buffer,999,area);
if (ok)
{
column = sscanf(buffer,"%d %d %d",&begin,&end,&sequence);
if (column > 1 && begin <= end)
{
begin--;
end--;
if (column==3)
sequence--;
else
sequence=0;
if (sequence && complement)
sequence++;
push_area(table,sequence,begin,end);
if (!sequence && complement)
push_area(table,1,complement -1 - end,complement -1 -begin);
}
if (column==1)
fprintf(stderr,"WARNING in mask file reading line %d\n",linecount);
}
} while (ok);
fprintf(stderr,"\nread %d masked areas from file\n",table->total);
strip_area_table(table);
fprintf(stderr,"strip to %d non overlaping areas\n",table->total);
return table;
}
char KMRK_isMasked(masked_area_table_t *mask,int32_t seq, int32_t position)
{
masked_area_t input;
int32_t result;
masked_area_list_t *list;
if (! mask || (seq >= mask->seqcount))
return 0;
list = mask->sequence[seq];
input.begin=position;
result = bsearch(&input,
list->area,
list->count,
sizeof(masked_area_t),
(int (*)(const void *, const void *))search_area) != NULL;
return result;
}
\ No newline at end of file
/*
* KMRK_mask.h
* repseek
*
* Created by Eric Coissac on 04/12/04.
* Copyright 2004 __MyCompanyName__. All rights reserved.
*
*/
#include <stdint.h>
#ifndef _KMRK_MASK_H_
#define _KMRK_MASK_H_
typedef struct {
int32_t begin;
int32_t end;
} masked_area_t;
typedef struct {
int32_t reserved;
int32_t count;
masked_area_t area[1];
} masked_area_list_t;
typedef struct {
int32_t seqcount;
int32_t total;
masked_area_list_t *sequence[1];
} masked_area_table_t;
masked_area_table_t *KMRK_ReadMaskArea(char *areafile,int32_t seqcount,int32_t complement);
char KMRK_isMasked(masked_area_table_t *mask,int32_t seq, int32_t position);
#endif
\ No newline at end of file
/**
* @file KMRK_merge_seeds.c
* @author Eric Coissac <coissac@inrialpes.fr>
* @date Wed Mar 3 11:15:57 2004
*
* @brief Merge function of overlapping seeds
*
*
*/
#include "KMRK_merge_seeds.h"
void KMRK_MergeSeeds(AllSeeds_type *AllSeeds,
int8_t opt_dir,
int8_t opt_inv)
{
int32_t i; /* the current seed */
int32_t j; /* the checked seed */
int32_t N; /* the kept ones */
Seed_type* seeds;
if(opt_dir){
seeds = AllSeeds->dirSeeds;
for(i=0, N=0 ;i<AllSeeds->nDirSeeds; i++){
if(seeds[i].pos1==-1) /* any seed at -1 is removed */
continue;
j=i+1;
while( (seeds[j].pos1!=-1) &&
(seeds[i].pos1!=-1) &&
(j < AllSeeds->nDirSeeds) &&
(seeds[j].pos1 < seeds[i].pos1+ seeds[i].length))
{
if(
((seeds[j].pos2 >= seeds[i].pos2) &&
(seeds[j].pos2 < seeds[i].pos2 + seeds[i].length)) || /* if the seeds are overlapping */
((seeds[j].pos2 + seeds[j].length >= seeds[i].pos2) &&
(seeds[j].pos2 + seeds[j].length < seeds[i].pos2 + seeds[i].length)))
{
if(seeds[j].length <= seeds[i].length)
seeds[j].pos1=seeds[j].pos2=seeds[j].length=-1; /* removed the smallest */
else
seeds[i].pos1=seeds[i].pos2=seeds[i].length=-1;
}
j++;
}
if(seeds[i].pos1 !=-1)
{ /* if this seed is not out, store it */
seeds[N].pos1 = seeds[i].pos1;
seeds[N].pos2 = seeds[i].pos2;
seeds[N].length = seeds[i].length;
N++;
}
}
AllSeeds->nFilteredDirSeeds += AllSeeds->nDirSeeds-N;
AllSeeds->nDirSeeds=N;
}
if(opt_inv){
seeds = AllSeeds->invSeeds;
for(i=0, N=0 ;i<AllSeeds->nInvSeeds; i++){
if(seeds[i].pos1==-1)
continue;
j=i+1;
while( (seeds[j].pos1!=-1 ) &&
(seeds[i].pos1!=-1 ) &&
(j < AllSeeds->nInvSeeds) &&
(seeds[j].pos1 < seeds[i].pos1+seeds[i].length))
{
if(
((seeds[j].pos2 >= seeds[i].pos2) && /* if the seeds are overlapping */
(seeds[j].pos2 < seeds[i].pos2+seeds[i].length)) ||
((seeds[j].pos2 + seeds[j].length >= seeds[i].pos2) &&
(seeds[j].pos2 + seeds[j].length < seeds[i].pos2+seeds[i].length)))
{
if(seeds[j].length <= seeds[i].length)
seeds[j].pos1=seeds[j].pos2=seeds[j].length=-1; /* removed the smallest */
else
seeds[i].pos1=seeds[i].pos2=seeds[i].length=-1;
}
j++;
}
if(seeds[i].pos1!=-1)
{ /* if this seed is not out, store it */
seeds[N].pos1 = seeds[i].pos1;
seeds[N].pos2 = seeds[i].pos2;
seeds[N].length = seeds[i].length;
N++;
}
}
AllSeeds->nFilteredInvSeeds += AllSeeds->nInvSeeds-N;
AllSeeds->nInvSeeds=N;
}
KMRK_compactSeeds(AllSeeds);
}
#ifndef KMRK_merge_seeds_h
#define KMRK_merge_seeds_h
#include "KMRK_Seeds.h"
void KMRK_MergeSeeds(AllSeeds_type *AllSeeds,
int8_t opt_dir,
int8_t opt_inv);
#endif /* KMRK_MergeSeeds */
SOURCES = KMRK.c \
KMRK_mask.c \
KMRK_merge_seeds.c \
KMRK_Seeds.c
SRCS=$(SOURCES)
OBJECTS= $(patsubst %.c,%.o,$(SOURCES))
LIBFILE= libKMRK.a
RANLIB= ranlib
include ../global.mk
all: $(LIBFILE)
clean:
rm -rf $(OBJECTS) $(LIBFILE)
$(LIBFILE): $(OBJECTS)
ar -cr $@ $?
$(RANLIB) $@
/******
file : memory.c
function : All about memory of the KMR, Seeds and Repeats
All MyMalloc() series is about follwoing how mauch memory has been used
created : 19 Sep 2003
modif : Oct 2003, Feb 2004
modif : Dec 2004 <EC> ; Corrected Memory declaration
author : amikezor
*****/
#include <stdio.h>
#include <stdlib.h>
#include "repseek_types.h"
#include "memory.h"
MemUsage Memory;
/*
Functions to count the memory usage all along
dybamic allocation and free
*/
void PrintMem(char *Comment){
extern MemUsage Memory;
fprintf(stderr,"\n%s\nMemory Usage\n\t* Max is: %d bytes, %.2f Kb, %.2f Mb\n\t* Cur is: %d bytes, %.2f Kb, %.2f Mb\n",
Comment,
Memory.Max, (float)Memory.Max/1024, (float)Memory.Max/(1024*1024),
Memory.Current, (float)Memory.Current/1024, (float)Memory.Current/(1024*1024));
}
void PrintMaxMem( void ){
extern MemUsage Memory;
if(Memory.Max < 1024)
fprintf(stderr,"Max Memory Usage.. %d bytes\n", Memory.Max);
else if(Memory.Max < 1024*1024)
fprintf(stderr,"Max Memory Usage.. %.2f Kilobytes\n", (float)Memory.Max/1024);
else if(Memory.Max < 1024*1024*1024)
fprintf(stderr,"Max Memory Usage.. %.2f Megabytes\n", (float)Memory.Max/(1024*1024));
else
fprintf(stderr,"Max Memory Usage.. %.2f Gigabytes\n", (float)Memory.Max/(1024*1024*1024));
}
void Update_Mem(int32_t Value){
extern MemUsage Memory;
Memory.Current += Value;
Memory.Max = <