fermi-lat / Likelihood

BSD 3-Clause "New" or "Revised" License
2 stars 1 forks source link

gtexpmap failure #70

Closed eacharles closed 4 years ago

eacharles commented 4 years ago

https://jira.slac.stanford.edu/projects/LK/issues/LK-133

eacharles commented 4 years ago

Ok. This was due a numerical precision issue in defining the ranges of integration in computing the Psf integral. The code tried to find upper and lower integration ranges by taking the cosine of the difference of two angles. When that value was small, the numerical precision caused it to try to integrate from 1 to 1 and caused this error. The basic issue is in the code block below.
I've replaced the test statement to use a full integration when roi_radius is very close to psi.


   double mup(std::cos(roi_radius + psi));
   double mum(std::cos(roi_radius - psi));
-   if (psi < roi_radius) {
+   if ( mum < 0.99 ) { 

      PsfIntegrand1 psfIntegrand1(sigma, gamma);
      firstIntegral = 
         st_facilities::GaussianQuadrature::dgaus8(psfIntegrand1, mum,
                                                   one, err, ierr);
   }
donhorner commented 4 years ago

Jean Ballet reports (using the Fermitools beta 1.3.0): "the gtexpmap call still fails with a very similar message Caught N13st_facilities18GaussianQuadrature15dgaus8ExceptionE at the top level: dgaus8 --- a and b are too nearly equal to allow normal integration. ans is set to 0 and ierr is set to -1."

eacharles commented 4 years ago

It is possible that there is still a case where this fails, but I fixed the previous failure.
Can we make sure that the release include the code that fixes this. (it is described in the GitHub issue).

-e

On May 29, 2020, at 10:28 AM, Don Horner notifications@github.com wrote:

Jean Ballet reports (using the Fermitools beta 1.3.0): "the gtexpmap call still fails with a very similar message Caught N13st_facilities18GaussianQuadrature15dgaus8ExceptionE at the top level: dgaus8 --- a and b are too nearly equal to allow normal integration. ans is set to 0 and ierr is set to -1."

— You are receiving this because you were assigned. Reply to this email directly, view it on GitHub https://github.com/fermi-lat/Likelihood/issues/70#issuecomment-636092770, or unsubscribe https://github.com/notifications/unsubscribe-auth/ADRIGIWVKUYWWDB37YJNRX3RT7WE7ANCNFSM4LUNLONA.

jasercion commented 4 years ago

I think I've identified the one remaining case where this would apply.

irfs/irfInterface/src/IPsf.cxx:281

I'll make the update and incorporate it into the rebuilt beta

eacharles commented 4 years ago

Nice job.

-e

On May 29, 2020, at 1:21 PM, jasercion notifications@github.com wrote:

I think I've identified the one remaining case where this would apply.

irfs/irfInterface/src/IPsf.cxx:181

I'll make the update and incorporate it into the rebuilt beta

— You are receiving this because you were assigned. Reply to this email directly, view it on GitHub https://github.com/fermi-lat/Likelihood/issues/70#issuecomment-636174567, or unsubscribe https://github.com/notifications/unsubscribe-auth/ADRIGIVB2QOAUKMHACN6HS3RUAKLRANCNFSM4LUNLONA.

jasercion commented 4 years ago

This PR was also hanging. It has now been merged. https://github.com/fermi-lat/irfs/pull/11

jballet commented 4 years ago

I have tested FT 1.3.1 and it does not fail any more. Of course this was a rare behavior so this single test cannot guarantee that gtexpmap does not fail in other rare cases. But it definitely solves this particular failure.

jasercion commented 4 years ago

@jballet Great, I'm going to close this for now then. If it happens again please reopen this issue.