sib-swiss / pftools3

A suite of tools to build and search generalized profiles
GNU General Public License v2.0
10 stars 7 forks source link

pfscale 2.3rev5d failure with certain profiles #31

Open beatrice79 opened 2 months ago

beatrice79 commented 2 months ago

The issue regards the writing of the R1 and R2 values. The part of the code concerns with the problem is wrprf.f lines 262-26 .

The problem manifests itself with the following command:

pfscale test_score.txt test_prf.txt N=21210388 P=0.000139130 Q=0.000 00047

With the output reporting:

....
At line 267 of file wrprf.f
Fortran runtime error: End of record
...

[test_score.txt] [test_prf.txt]

In this peculiar case, the problem concerns the writing of the R1 value (RNOP(I2,I1) in the code) which is 17.6187992 (If printed with 7 digits to the right of the decimal point)

Rewriting the IF/ELSE clause of wrprf.f lines 262-26 in pseudo code:

IF ((|R1| ≤ 10 AND |R1| > 0.0001) OR (R1 = 0)):
Write(CHRP(I2),'(F10.7)') R1
ELSE:
Write(CHRP(I2),*) R1

A few comments/interrogations:

smoretti commented 2 months ago

Found the source code of version 2.3.2 if you want it (also earlier ones).

But the differences between wrprf.f versions do not look striking:

2c2
< * $Id: wrprf.f,v 2.8 2003/07/03 13:08:58 vflegel Exp $
---
> * $Id: wrprf.f,v 2.9 2003/11/18 10:52:29 vflegel Exp $
207c207
<          JP=Lblnk(CPAR)
---
>          JP=Lblnk(CPAR)+1
212c212
<       JB=JB+JP+1
---
>       JB=JB+JP
269,270c269,270
<             Write(CPAR(JP+1:),*)'R',I2,'=',CHRP(I2),';'
<             JP=Lblnk(CPAR)
---
>             Write(CPAR(JP+1:),*)'R',I2,'=',CHRP(I2),'; '
>             JP=Lblnk(CPAR)+1
275c275
<          JB=JB+JP+1
---
>          JB=JB+JP

Of course, like you mentioned, no information about this 10 upper limit.

Did you try to replace 10 by, for example 20, to see if the result makes sense?

beatrice79 commented 1 month ago

Yes I did try replacing 10 by 20 and the profile gets produced is fine, it has an R1 value of 17.6187992 when the one produced by pfscale 2.3 rev 2 has an R1 value of 17.6187973 - not much of a difference.

smoretti commented 1 month ago

Did you see some differences with pfscale 2.3 rev 2 I missed?

beatrice79 commented 1 month ago

Did you see some differences with pfscale 2.3 rev 2 I missed? I don't know about relevant differences. pfscale.f makes use of:

Files used by pfscale.f (even transitively) that differ between 2.3 rev 2 and rev 5.d are: