Papara is "bailing out"

[fangly]

Hi,

I am using Papara 2.4, trying to align a few sequences into Greengenes. I get a "bailing out" error / core dump:

papara called as:
papara -t gg_97_otus_unannotated_indented.tree -s gg_97_otus_aln.phylip -q all_seqs_clean.fna 
gap rate: 20605051 147592492
gap rate: 0.139608
rate matrix: [2,2]((-0.139608,0.139608),(0.860392,-0.860392))
p: [2,2]((0.986715,0.0132854),(0.0818772,0.918123))
meeeeep: 0 0
terminate called after throwing an instance of 'std::runtime_error'
  what():  bailing out.
Aborted (core dumped)

This seems to refer to the function gap_posterior, in which both v1 and v2 are 0, which means that there would be a division/0 problem when calculating v:

    static inline double gap_posterior( double v1, double v2 ) {
        assert( pgap_model.is_valid_ptr() );
        //return v1 / (v1 + v2);

        v1 *= 1 - pgap_model->gap_freq();
    //  throw std::runtime_error( "i think there is an error in this function. why v1 in the next line?");
        v2 *= pgap_model->gap_freq();

        float v = float(v1 / (v1 + v2));

        if( v != v ) {
                std::cerr << "meeeeep: " << v1 << " " << v2 << "\n";

                throw std::runtime_error("bailing out.");
        }

        return v;
    }

Any idea what is going on? I could provide some data if necessary.

[fangly]

This patch is working for me, but I am not 100% sure of its correctness:

/srv/sw/papara/2.4$ diff pvec.h.ori pvec.h
549c549,551
<       float v = float(v1 / (v1 + v2));
---
>   float v;
>   if (v1 != 0) {
>       v = float(v1 / (v1 + v2));
551,552c553,554
<       if( v != v ) {
<           std::cerr << "meeeeep: " << v1 << " " << v2 << "\n";
---
>           if( v != v ) {
>               std::cerr << "meeeeep: " << v1 << " " << v2 << "\n";
554,555c556,560
<           throw std::runtime_error("bailing out.");
<       }
---
>           throw std::runtime_error("bailing out.");
>           }
>   } else {
>       v = v1;
>   }