From 32e4df4d446ac72939e17e5553d715e8c1e80efc Mon Sep 17 00:00:00 2001 From: Kevin Moran <39282665+Thernn88@users.noreply.github.com> Date: Fri, 4 Sep 2026 19:46:31 -0500 Subject: [PATCH 1/2] Add lazy bathalign, and --nofs to skip frameshift calibration BATH dropped HMMER's hmmalign driver but kept the whole tracealign engine (tracealign.c, p7_tracealign_*) compiled into libhmmer.a. bathalign is a 288-line driver over that existing engine -- no new alignment code. Alphabet-generic, as the engine is; BATH3/f is a superset of HMMER3/f, so it reads both .bhmm and .hmm. Amino output is byte-identical to HMMER 3.4 hmmalign across --mapali, --outformat, --trim, and Stockholm and afa output, for models from bathbuild and hmmbuild alike. Note that bathbuild itself is amino-only, so non-amino models must come from hmmbuild. c2ba68e made frameshift tau calculation unconditional and dropped --fs, but left the gating in place: evalues.c still guards FS3/FS5 on cfg_b->fs and p7_hmmfile.c still guards the header lines on hmm->fs. --nofs re-exposes that switch against the new default. Alignment never reads the frameshift stats, so for a model built to be aligned against, the calibration is pure waste -- and it is more than half of bathbuild's runtime. Where the model exists only for alignment, building it dominates the cost of aligning. 400 dipteraODB12 alignments, serial, 3 runs each, in seconds: full 43.83 / 43.32 / 43.38 nofs 15.26 / 15.50 / 15.57 A --nofs model is byte-identical to a full one apart from the omitted FS3/FS5 STATS and FRAMESHIFT PROB lines, and gives identical bathsearch results when --fs/--fsonly are not used. When they are, bathsearch fails cleanly on the unset stats rather than mis-scoring. --- src/Makefile.in | 2 + src/bathalign.c | 288 ++++++++++++++++++++++++++++++++++++++++++++++++ src/bathbuild.c | 6 + 3 files changed, 296 insertions(+) create mode 100644 src/bathalign.c diff --git a/src/Makefile.in b/src/Makefile.in index 6beceb19..152ecef5 100644 --- a/src/Makefile.in +++ b/src/Makefile.in @@ -49,12 +49,14 @@ MYINCDIRS = -I../${ESLDIR} \ PROGS = bathsearch\ + bathalign\ bathbuild\ bathconvert\ bathfetch\ bathstat PROGOBJS = bathsearch.o\ + bathalign.o\ bathbuild.o\ bathconvert.o\ bathfetch.o\ diff --git a/src/bathalign.c b/src/bathalign.c new file mode 100644 index 00000000..54fb8278 --- /dev/null +++ b/src/bathalign.c @@ -0,0 +1,288 @@ +/* bathalign: align protein sequences to a BATH profile HMM + */ +#include "p7_config.h" + +#include +#include +#include + +#include "easel.h" +#include "esl_alphabet.h" +#include "esl_getopts.h" +#include "esl_msa.h" +#include "esl_msafile.h" +#include "esl_sq.h" +#include "esl_sqio.h" +#include "esl_vectorops.h" + +#include "hmmer.h" + +static int map_alignment(const char *msafile, const P7_HMM *hmm, ESL_SQ ***ret_sq, P7_TRACE ***ret_tr, int *ret_ntot); + + +#define ALPHOPTS "--amino,--dna,--rna" /* Exclusive options for alphabet choice */ + +static ESL_OPTIONS options[] = { + /* name type default env range toggles reqs incomp help docgroup*/ + { "-h", eslARG_NONE, FALSE, NULL, NULL, NULL, NULL, NULL, "show brief help on version and usage", 1 }, + { "-o", eslARG_OUTFILE, NULL, NULL, NULL, NULL, NULL, NULL, "output alignment to file , not stdout", 1 }, + + { "--mapali", eslARG_INFILE, NULL, NULL, NULL, NULL, NULL, NULL, "include alignment in file (same ali that HMM came from)", 2 }, + { "--trim", eslARG_NONE, FALSE, NULL, NULL, NULL, NULL, NULL, "trim terminal tails of nonaligned residues from alignment", 2 }, + { "--amino", eslARG_NONE, FALSE, NULL, NULL, ALPHOPTS, NULL, NULL, "assert , both protein: no autodetection", 2 }, + { "--dna", eslARG_NONE, FALSE, NULL, NULL, ALPHOPTS, NULL, NULL, "assert , both DNA: no autodetection", 2 }, + { "--rna", eslARG_NONE, FALSE, NULL, NULL, ALPHOPTS, NULL, NULL, "assert , both RNA: no autodetection", 2 }, + { "--informat", eslARG_STRING, NULL, NULL, NULL, NULL, NULL, NULL, "assert is in format : no autodetection", 2 }, + { "--outformat", eslARG_STRING, "Stockholm", NULL, NULL, NULL, NULL, NULL, "output alignment in format ", 2 }, + { 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 }, +}; + +static char usage[] = "[-options] "; +static char banner[] = "align protein sequences to a profile HMM"; + + + + +static void +cmdline_failure(char *argv0, char *format, ...) +{ + va_list argp; + printf("\nERROR: "); + va_start(argp, format); + vfprintf(stdout, format, argp); + va_end(argp); + esl_usage(stdout, argv0, usage); + printf("\nTo see more help on available options, do %s -h\n\n", argv0); + exit(1); +} + +static void +cmdline_help(char *argv0, ESL_GETOPTS *go) +{ + p7_banner (stdout, argv0, banner); + esl_usage (stdout, argv0, usage); + puts("\nBasic options:"); + esl_opt_DisplayHelp(stdout, go, 1, 2, 80); + puts("\nLess common options:"); + esl_opt_DisplayHelp(stdout, go, 2, 2, 80); + puts("\nSequence input formats include: FASTA, EMBL, GenBank, UniProt"); + puts("Alignment output formats include: Stockholm, Pfam, A2M, PSIBLAST\n"); + exit(0); +} + + + +int +main(int argc, char **argv) +{ + ESL_GETOPTS *go = NULL; /* application configuration */ + char *hmmfile = NULL; /* HMM file name */ + char *seqfile = NULL; /* sequence file name */ + char *mapfile = NULL; /* optional mapped MSA file name */ + int infmt = eslSQFILE_UNKNOWN; + int outfmt = eslMSAFILE_STOCKHOLM; + P7_HMMFILE *hfp = NULL; /* open HMM file */ + ESL_SQFILE *sqfp = NULL; /* open sequence file */ + char *outfile = NULL; /* output filename */ + FILE *ofp = stdout; /* output stream */ + ESL_SQ **sq = NULL; /* array of sequences */ + void *p = NULL; /* tmp ptr for reallocation */ + int nseq = 0; /* # of sequences in */ + int mapseq = 0; /* # of sequences in mapped MSA */ + int totseq = 0; /* # of seqs in all sources */ + ESL_ALPHABET *abc = NULL; /* alphabet (set from the HMM file)*/ + P7_HMM *hmm = NULL; + P7_TRACE **tr = NULL; /* array of tracebacks */ + ESL_MSA *msa = NULL; /* resulting multiple alignment */ + int msaopts = 0; /* flags to p7_tracealign_Seqs() */ + int idx; /* counter over seqs, traces */ + int status; /* easel/hmmer return code */ + char errbuf[eslERRBUFSIZE]; + + /* Parse the command line + */ + go = esl_getopts_Create(options); + if (esl_opt_ProcessCmdline(go, argc, argv) != eslOK) cmdline_failure(argv[0], "Failed to parse command line: %s\n", go->errbuf); + if (esl_opt_VerifyConfig(go) != eslOK) cmdline_failure(argv[0], "Error in configuration: %s\n", go->errbuf); + if (esl_opt_GetBoolean(go, "-h") ) cmdline_help (argv[0], go); + if (esl_opt_ArgNumber(go) != 2) cmdline_failure(argv[0], "Incorrect number of command line arguments.\n"); + + hmmfile = esl_opt_GetArg(go, 1); + seqfile = esl_opt_GetArg(go, 2); + + if (strcmp(hmmfile, "-") == 0 && strcmp(seqfile, "-") == 0) + cmdline_failure(argv[0], "Either or may be '-' (to read from stdin), but not both.\n"); + + msaopts |= p7_ALL_CONSENSUS_COLS; /* default as of 3.1 */ + if (esl_opt_GetBoolean(go, "--trim")) msaopts |= p7_TRIM; + + /* If caller declared an input format, decode it + */ + if (esl_opt_IsOn(go, "--informat")) { + infmt = esl_sqio_EncodeFormat(esl_opt_GetString(go, "--informat")); + if (infmt == eslSQFILE_UNKNOWN) cmdline_failure(argv[0], "%s is not a recognized input sequence file format\n", esl_opt_GetString(go, "--informat")); + } + + /* Determine output alignment file format */ + outfmt = esl_msafile_EncodeFormat(esl_opt_GetString(go, "--outformat")); + if (outfmt == eslMSAFILE_UNKNOWN) cmdline_failure(argv[0], "%s is not a recognized output MSA file format\n", esl_opt_GetString(go, "--outformat")); + + /* Open output stream */ + if ( (outfile = esl_opt_GetString(go, "-o")) != NULL) + { + if ((ofp = fopen(outfile, "w")) == NULL) + cmdline_failure(argv[0], "failed to open -o output file %s for writing\n", outfile); + } + + + /* If caller forced an alphabet on us, create the one the caller wants + */ + if (esl_opt_GetBoolean(go, "--amino")) abc = esl_alphabet_Create(eslAMINO); + else if (esl_opt_GetBoolean(go, "--dna")) abc = esl_alphabet_Create(eslDNA); + else if (esl_opt_GetBoolean(go, "--rna")) abc = esl_alphabet_Create(eslRNA); + + /* Read one HMM, and make sure there's only one. + */ + status = p7_hmmfile_OpenE(hmmfile, NULL, &hfp, errbuf); + if (status == eslENOTFOUND) p7_Fail("File existence/permissions problem in trying to open HMM file %s.\n%s\n", hmmfile, errbuf); + else if (status == eslEFORMAT) p7_Fail("File format problem in trying to open HMM file %s.\n%s\n", hmmfile, errbuf); + else if (status != eslOK) p7_Fail("Unexpected error %d in opening HMM file %s.\n%s\n", status, hmmfile, errbuf); + + status = p7_hmmfile_Read(hfp, &abc, &hmm); + if (status == eslEFORMAT) p7_Fail("Bad file format in HMM file %s:\n%s\n", hfp->fname, hfp->errbuf); + else if (status == eslEINCOMPAT) p7_Fail("HMM in %s is not in the expected %s alphabet\n", hfp->fname, esl_abc_DecodeType(abc->type)); + else if (status == eslEOF) p7_Fail("Empty HMM file %s? No HMM data found.\n", hfp->fname); + else if (status != eslOK) p7_Fail("Unexpected error in reading HMMs from %s\n", hfp->fname); + + status = p7_hmmfile_Read(hfp, &abc, NULL); + if (status != eslEOF) p7_Fail("HMM file %s does not contain just one HMM\n", hfp->fname); + p7_hmmfile_Close(hfp); + + + /* We're going to build up two arrays: sequences and traces. + * If --mapali option is chosen, the first set of sequences/traces is from the provided alignment + */ + if ( (mapfile = esl_opt_GetString(go, "--mapali")) != NULL) + { + map_alignment(mapfile, hmm, &sq, &tr, &mapseq); + } + totseq = mapseq; + + /* Read digital sequences into an array (possibly concat'ed onto mapped seqs) + */ + status = esl_sqfile_OpenDigital(abc, seqfile, infmt, NULL, &sqfp); + if (status == eslENOTFOUND) p7_Fail("Failed to open sequence file %s for reading\n", seqfile); + else if (status == eslEFORMAT) p7_Fail("Sequence file %s is empty or misformatted\n", seqfile); + else if (status != eslOK) p7_Fail("Unexpected error %d opening sequence file %s\n", status, seqfile); + + ESL_RALLOC(sq, p, sizeof(ESL_SQ *) * (totseq + 1)); + sq[totseq] = esl_sq_CreateDigital(abc); + nseq = 0; + while ((status = esl_sqio_Read(sqfp, sq[totseq+nseq])) == eslOK) + { + nseq++; + ESL_RALLOC(sq, p, sizeof(ESL_SQ *) * (totseq+nseq+1)); + sq[totseq+nseq] = esl_sq_CreateDigital(abc); + } + if (status == eslEFORMAT) esl_fatal("Parse failed (sequence file %s):\n%s\n", + sqfp->filename, esl_sqfile_GetErrorBuf(sqfp)); + else if (status != eslEOF) esl_fatal("Unexpected error %d reading sequence file %s", status, sqfp->filename); + esl_sqfile_Close(sqfp); + totseq += nseq; + + + /* Remaining initializations, including trace array allocation + */ + ESL_RALLOC(tr, p, sizeof(P7_TRACE *) * totseq); + for (idx = mapseq; idx < totseq; idx++) + tr[idx] = p7_trace_CreateWithPP(); + + p7_tracealign_computeTraces(hmm, sq, mapseq, totseq - mapseq, tr); + + p7_tracealign_Seqs(sq, tr, totseq, hmm->M, msaopts, hmm, &msa); + + esl_msafile_Write(ofp, msa, outfmt); + + for (idx = 0; idx <= totseq; idx++) esl_sq_Destroy(sq[idx]); /* including sq[nseq] because we overallocated */ + for (idx = 0; idx < totseq; idx++) p7_trace_Destroy(tr[idx]); + free(sq); + free(tr); + esl_msa_Destroy(msa); + p7_hmm_Destroy(hmm); + if (ofp != stdout) fclose(ofp); + esl_alphabet_Destroy(abc); + esl_getopts_Destroy(go); + return eslOK; + + ERROR: + return status; +} + + + +/***************************************************************** + * Internal functions used by main and API + *****************************************************************/ + +static int +map_alignment(const char *msafile, const P7_HMM *hmm, ESL_SQ ***ret_sq, P7_TRACE ***ret_tr, int *ret_ntot) +{ + ESL_SQ **sq = NULL; + P7_TRACE **tr = NULL; + ESL_MSAFILE *afp = NULL; + ESL_MSA *msa = NULL; + ESL_ALPHABET *abc = (ESL_ALPHABET *) hmm->abc; /* removing const'ness to make compiler happy. Safe. */ + int *matassign = NULL; + uint32_t chksum = 0; + int i,k; + int status; + + status = esl_msafile_Open(&abc, msafile, NULL, eslMSAFILE_UNKNOWN, NULL, &afp); + if (status != eslOK) esl_msafile_OpenFailure(afp, status); + + status = esl_msafile_Read(afp, &msa); + if (status != eslOK) esl_msafile_ReadFailure(afp, status); + + if (! (hmm->flags & p7H_CHKSUM) ) esl_fatal("HMM has no checksum. --mapali unreliable without it."); + if (! (hmm->flags & p7H_MAP) ) esl_fatal("HMM has no map. --mapali can't work without it."); + esl_msa_Checksum(msa, &chksum); + if (hmm->checksum != chksum) esl_fatal("--mapali MSA %s isn't same as the one HMM came from (checksum mismatch)", msafile); + + ESL_ALLOC(sq, sizeof(ESL_SQ *) * msa->nseq); + ESL_ALLOC(tr, sizeof(P7_TRACE *) * msa->nseq); + ESL_ALLOC(matassign, sizeof(int) * (msa->alen + 1)); + + esl_vec_ISet(matassign, msa->alen+1, 0); + for (k = 1; k <= hmm->M; k++) matassign[hmm->map[k]] = 1; + + p7_trace_FauxFromMSA(msa, matassign, p7_DEFAULT, tr); + + /* The 'faux' core traces constructed by FauxFromMSA() may contain + * D->I and I->D transitions. They may *only* now be passed to + * p7_tracealign_Seqs(), which can deal with these 'illegal' + * transitions, in order to exactly reproduce the input --mapali + * alignment. + */ + + for (i = 0; i < msa->nseq; i++) + esl_sq_FetchFromMSA(msa, i, &(sq[i])); + + *ret_ntot = msa->nseq; + *ret_tr = tr; + *ret_sq = sq; + + esl_msafile_Close(afp); + esl_msa_Destroy(msa); + free(matassign); + return eslOK; + + ERROR: + *ret_ntot = 0; + *ret_tr = NULL; + *ret_sq = NULL; + if (afp != NULL) esl_msafile_Close(afp); + if (msa != NULL) esl_msa_Destroy(msa); + if (matassign != NULL) free(matassign); + return status; +} + diff --git a/src/bathbuild.c b/src/bathbuild.c index 7aa5ae5f..2b1ef3f7 100644 --- a/src/bathbuild.c +++ b/src/bathbuild.c @@ -67,6 +67,7 @@ static ESL_OPTIONS options[] = { { "-o", eslARG_OUTFILE,FALSE, NULL, NULL, NULL, NULL, NULL, "direct summary output to file , not stdout", 1 }, { "-O", eslARG_OUTFILE,FALSE, NULL, NULL, NULL, NULL, NULL, "resave annotated, possibly modified MSA to file ", 1 }, { "--ct", eslARG_INT, "1", NULL, NULL, NULL, NULL, NULL, "use alt genetic code of NCBI transl table ", 1 }, + { "--nofs", eslARG_NONE, FALSE, NULL, NULL, NULL, NULL, NULL, "skip frameshift stats; model is align-only, not searchable", 1 }, /* Alternate model construction strategies */ { "--fast", eslARG_NONE, "default",NULL, NULL, CONOPTS, NULL, NULL, "assign cols w/ >= symfrac residues as consensus", 3 }, @@ -264,6 +265,7 @@ output_header(const ESL_GETOPTS *go, const struct cfg_s *cfg) if (fprintf(cfg->ofp, "# input file: %s\n", cfg->infile) < 0) ESL_EXCEPTION_SYS(eslEWRITE, "write failed"); if (fprintf(cfg->ofp, "# output HMM file: %s\n", cfg->hmmfile) < 0) ESL_EXCEPTION_SYS(eslEWRITE, "write failed"); + if (esl_opt_GetBoolean(go, "--nofs") && fprintf(cfg->ofp, "# frameshift stats calculated: NO (align-only model)\n") < 0) ESL_EXCEPTION_SYS(eslEWRITE, "write failed"); if (esl_opt_IsUsed(go, "-n") && fprintf(cfg->ofp, "# name (the single) HMM: %s\n", esl_opt_GetString(go, "-n")) < 0) ESL_EXCEPTION_SYS(eslEWRITE, "write failed"); if (esl_opt_IsUsed(go, "-o") && fprintf(cfg->ofp, "# output directed to file: %s\n", esl_opt_GetString(go, "-o")) < 0) ESL_EXCEPTION_SYS(eslEWRITE, "write failed"); @@ -608,6 +610,10 @@ usual_master(const ESL_GETOPTS *go, struct cfg_s *cfg) /* special arguments for hmmbuild */ info[i].bld->fsprob = p7P_FSPROB; + /* --nofs: skip the FS3/FS5 tau calibrations. Cheaper, but the resulting + * model carries no frameshift stats, so bathsearch will reject it. + * Intended for tools that only need the core profile (e.g. bathalign). */ + info[i].bld->fs = (! esl_opt_GetBoolean(go, "--nofs")); info[i].bld->w_len = (go != NULL && esl_opt_IsOn (go, "--w_length")) ? esl_opt_GetInteger(go, "--w_length"): -1; info[i].bld->w_beta = (go != NULL && esl_opt_IsOn (go, "--w_beta")) ? esl_opt_GetReal (go, "--w_beta") : p7_DEFAULT_WINDOW_BETA; if ( info[i].bld->w_beta < 0 || info[i].bld->w_beta > 1 ) esl_fatal("Invalid window-length beta value\n"); From 1ae94c1f5a8a0268306593590a323dfae38c7722 Mon Sep 17 00:00:00 2001 From: Kevin Moran <39282665+Thernn88@users.noreply.github.com> Date: Thu, 17 Sep 2026 14:30:38 -0500 Subject: [PATCH 2/2] bathbuild: make frameshift calibration opt-in via --fs --fs is back, opt-in, default off. Partial revert of c2ba68e7: same option name and default as before April, keeping that commit's bathconvert/bathfetch work, which lets a model get FS3/FS5 added later instead of rebuilt. Help text and summary line fixed too. This changes the default on main. Anyone doing bathbuild then bathsearch --fs will need to add --fs or run bathconvert. i9-optional-annotation.pl built a model without --fs and then searched it with --fs; it now builds with --fs. bathsearch had three --fs checks all saying "not formated for this version bathsearch". Under the new default that error usually just means the model was built without --fs, so I merged them into one check that says what's missing and how to fix it. Co-Authored-By: Claude Opus 5 --- src/bathbuild.c | 10 ++++------ src/bathsearch.c | 7 ++++--- testsuite/i9-optional-annotation.pl | 2 +- 3 files changed, 9 insertions(+), 10 deletions(-) diff --git a/src/bathbuild.c b/src/bathbuild.c index 2b1ef3f7..2228eef0 100644 --- a/src/bathbuild.c +++ b/src/bathbuild.c @@ -67,7 +67,7 @@ static ESL_OPTIONS options[] = { { "-o", eslARG_OUTFILE,FALSE, NULL, NULL, NULL, NULL, NULL, "direct summary output to file , not stdout", 1 }, { "-O", eslARG_OUTFILE,FALSE, NULL, NULL, NULL, NULL, NULL, "resave annotated, possibly modified MSA to file ", 1 }, { "--ct", eslARG_INT, "1", NULL, NULL, NULL, NULL, NULL, "use alt genetic code of NCBI transl table ", 1 }, - { "--nofs", eslARG_NONE, FALSE, NULL, NULL, NULL, NULL, NULL, "skip frameshift stats; model is align-only, not searchable", 1 }, + { "--fs", eslARG_NONE, FALSE, NULL, NULL, NULL, NULL, NULL, "calculate FS3/FS5 stats, for bathsearch --fs/--fsonly", 1 }, /* Alternate model construction strategies */ { "--fast", eslARG_NONE, "default",NULL, NULL, CONOPTS, NULL, NULL, "assign cols w/ >= symfrac residues as consensus", 3 }, @@ -265,7 +265,7 @@ output_header(const ESL_GETOPTS *go, const struct cfg_s *cfg) if (fprintf(cfg->ofp, "# input file: %s\n", cfg->infile) < 0) ESL_EXCEPTION_SYS(eslEWRITE, "write failed"); if (fprintf(cfg->ofp, "# output HMM file: %s\n", cfg->hmmfile) < 0) ESL_EXCEPTION_SYS(eslEWRITE, "write failed"); - if (esl_opt_GetBoolean(go, "--nofs") && fprintf(cfg->ofp, "# frameshift stats calculated: NO (align-only model)\n") < 0) ESL_EXCEPTION_SYS(eslEWRITE, "write failed"); + if (fprintf(cfg->ofp, "# frameshift stats calculated: %s\n", (esl_opt_GetBoolean(go, "--fs") ? "YES" : "NO")) < 0) ESL_EXCEPTION_SYS(eslEWRITE, "write failed"); if (esl_opt_IsUsed(go, "-n") && fprintf(cfg->ofp, "# name (the single) HMM: %s\n", esl_opt_GetString(go, "-n")) < 0) ESL_EXCEPTION_SYS(eslEWRITE, "write failed"); if (esl_opt_IsUsed(go, "-o") && fprintf(cfg->ofp, "# output directed to file: %s\n", esl_opt_GetString(go, "-o")) < 0) ESL_EXCEPTION_SYS(eslEWRITE, "write failed"); @@ -610,10 +610,8 @@ usual_master(const ESL_GETOPTS *go, struct cfg_s *cfg) /* special arguments for hmmbuild */ info[i].bld->fsprob = p7P_FSPROB; - /* --nofs: skip the FS3/FS5 tau calibrations. Cheaper, but the resulting - * model carries no frameshift stats, so bathsearch will reject it. - * Intended for tools that only need the core profile (e.g. bathalign). */ - info[i].bld->fs = (! esl_opt_GetBoolean(go, "--nofs")); + //frameshift stats are opt-in, as in bathsearch; bathconvert adds them later if needed + info[i].bld->fs = esl_opt_GetBoolean(go, "--fs"); info[i].bld->w_len = (go != NULL && esl_opt_IsOn (go, "--w_length")) ? esl_opt_GetInteger(go, "--w_length"): -1; info[i].bld->w_beta = (go != NULL && esl_opt_IsOn (go, "--w_beta")) ? esl_opt_GetReal (go, "--w_beta") : p7_DEFAULT_WINDOW_BETA; if ( info[i].bld->w_beta < 0 || info[i].bld->w_beta > 1 ) esl_fatal("Invalid window-length beta value\n"); diff --git a/src/bathsearch.c b/src/bathsearch.c index 12224f04..fd0cf53d 100644 --- a/src/bathsearch.c +++ b/src/bathsearch.c @@ -745,9 +745,10 @@ serial_master(ESL_GETOPTS *go, struct cfg_s *cfg) om = NULL; /* optimized query profile */ if(esl_opt_IsUsed(go, "--fs") || esl_opt_IsUsed(go, "--fsonly")) { //check that HMM is properly formated for bathsearch - if(!(hmm->fsprob && hmm->ct)) p7_Fail("HMM file %s not formated for this version bathsearch. Please run 'bathconvert new_file.bhmm old_file.bhmm'.\n", cfg->queryfile); - if( hmm->evparam[p7_FTAUFS3] == p7_EVPARAM_UNSET ) p7_Fail("HMM file %s not formated for this version bathsearch. Please run 'bathconvert new_file.bhmm old_file.bhmm'.\n", cfg->queryfile); - if( hmm->evparam[p7_FTAUFS5] == p7_EVPARAM_UNSET ) p7_Fail("HMM file %s not formated for this version bathsearch. Please run 'bathconvert new_file.bhmm old_file.bhmm'.\n", cfg->queryfile); + if( !(hmm->fsprob && hmm->ct) || + hmm->evparam[p7_FTAUFS3] == p7_EVPARAM_UNSET || + hmm->evparam[p7_FTAUFS5] == p7_EVPARAM_UNSET ) + p7_Fail("HMM file %s has no frameshift statistics, which --fs/--fsonly require.\nRebuild with 'bathbuild --fs', or add them with 'bathconvert new_file.bhmm %s'.\n", cfg->queryfile, cfg->queryfile); } else { diff --git a/testsuite/i9-optional-annotation.pl b/testsuite/i9-optional-annotation.pl index 46940aac..bb0f0366 100755 --- a/testsuite/i9-optional-annotation.pl +++ b/testsuite/i9-optional-annotation.pl @@ -58,7 +58,7 @@ BEGIN close ALI1; close SEQ1; -@output = `$builddir/src/bathbuild $tmppfx.bhmm $tmppfx.sto 2>&1`; +@output = `$builddir/src/bathbuild --fs $tmppfx.bhmm $tmppfx.sto 2>&1`; if ($? != 0) { die "FAIL: bathbuild failed\n"; } @output = `$builddir/src/bathsearch --fs --tblout $tmppfx.tbl $tmppfx.bhmm $tmppfx.seq 2>&1`;