Skip to content

Fix for #68#117

Merged
smaret merged 13 commits intolime-rt:masterfrom
allegroLeiden:fix_68
Jul 13, 2016
Merged

Fix for #68#117
smaret merged 13 commits intolime-rt:masterfrom
allegroLeiden:fix_68

Conversation

@allegroLeiden
Copy link
Copy Markdown
Contributor

@allegroLeiden allegroLeiden commented May 10, 2016

Fix for #68

In parseInput() the user-supplied magfield() function is tested and if the result is
non-zero, the polarization flag is set. (The default magfield function makes no changes
in the values of B, which are initialized to zero within parseInput().) Personally I
think it would be much simpler and certainly much more robust in SE terms just to have
the user set the flag directly. For the moment I have just made the test a little more
comprehensive, testing the magnitude of the B field returned rather than simply its X
component. The default magfield() function now also returns a value - again it is not
good SE practice not to return any value.
I changed the test value for (*img)[i].freq on line 148 of aux.c from 0 to the default of
-1.
This is just a cosmetic change, but 'angle' is another of the misleading names with
which Lime abounds.
- I changed stokesangles() to fix lime-rt#68. The image rotation matrix is now used to rotate
  the B field to the observer frame before the XY-plane angles are used to calculate the
  Stokes emissions, in a way also that uses no transcendental function calls. The changed
  interface was propagated out to sourceFunc_pol() and from there to traceray().

- I rewrote sourceFunc_pol() to be clearer and simpler.

- In traceray(), sourceFunc_pol() was being called for each of Stokes I, Q and U but it
  only needs to be called once for all. This is now fixed. I also make use now of
  calcSourceFn() in the polarisation block.
As pointed out by Sergey Parfenov, the expression for Stokes I was taken from a
paper which has been subsequently shown to be wrong by a factor of 2. This has now
been corrected here.
@smaret smaret added the bug label Jun 13, 2016
@smaret smaret added this to the Release 1.6 milestone Jun 13, 2016
...to reflect the corrected expression for sigma2 as prescribed by Ade et al., A&A 576, A105 (2015).
@tlunttil
Copy link
Copy Markdown
Contributor

tlunttil commented Jul 6, 2016

This PR looks good to me, but I have one suggestion/question:
Why is the check for polarization images done with an ugly and possibly error-prone test for a non-zero B field strength at a single point instead of (well, more like in addition to) par->polarization? The manual tells to set par->polarization if Stokes I/Q/U images are needed - there's nothing to indicate that setting B field alone is enough (and IMO, it shouldn't be). There's even a comment in the code saying that checking B!=0 isn't a good test, but is there some deep reason why it is still there?

A couple of very small cosmetic changes. I'm not sure if this PR is the right place to change these, or maybe these have already been done somewhere. They just caught my eye when I was reading the code.

  • There's a magic-number AU-to-parsec conversion on line 170 of aux.c - it would be cleaner if such things were done with constants defined in lime.h.
  • Function invSqrt (strarting line 338 of aux.c) isn't used (was it ever?) and should be removed. Anyway, the way it's done is not even a good idea on modern computers.

@allegroLeiden
Copy link
Copy Markdown
Contributor Author

This is why we need #77! I agree with you entirely @tlunttil about the test for B, which I have here only made marginally less evil than it was originally, but I had believed that par->polarization was an 'unsettable inpar' because whatever value the user might give it becomes overwritten at aux.c:141 or 155 (line nos in the PR version of aux.c). Apparently the doco believes otherwise.

I'm happy to change the code accordingly. I would (i) make par->polarization the sole means of controlling the polarization; (ii) still test the B field but only issue a warning if(par->polarization && normBSquared<=0).

I agree about the other 2 points but this is not the only place in lime with a hard-wired value. I had in mind at some point (not necessarily pre-1.6) going through systematically and fixing them all. If it makes you happy though I will fix them here. 😄

@tlunttil
Copy link
Copy Markdown
Contributor

tlunttil commented Jul 7, 2016

I think the cosmetic changes are better left to another PR (and yes, probably post-1.6). It's better to do a general code clean-up all at once instead of piecemeal in PRs that are meant for something else. I guess we'd better make new issue for that.

The documentation says that Stokes-parameters output for continuum images is switched on by setting par->polarization. The code did not
agree with this however: it tested the value returned by magfield() and if this was non-zero, it set par->polarization. Whatever value the
user provided was overwritten in all circumstances. This has now been changed. The user's choice now has the effect the documentation
describes. The test of magfield(0,0,0) is still performed, but only now issues a warning in the case of a zero value.
@tlunttil
Copy link
Copy Markdown
Contributor

I don't think you need to (and you shouldn't) remove setting the default values at the beginning of parseInput. In your new commit you don't seem to call anywhere input (defined in model.c), where all the parameters are usually set.

Also: if (normBSquared > 0.) on line 61??

@allegroLeiden
Copy link
Copy Markdown
Contributor Author

@tlunttil the goal posts have moved considerably since Reinhold Schaaf's PR was merged in. Pseudo-code is now

main.c:main(){
  main.c:initParImg(){
    // sets inpars defaults
    // reads them from input()
  }

  main.c:run(){
    aux.c:parseInput(){
      // does the remaining stuff with the pars
    }
    // does the remaining stuff which used to be in main()
  }
  // etc
}

Also: if (normBSquared > 0.) on line 61??

Sorry, what's the issue?

@allegroLeiden
Copy link
Copy Markdown
Contributor Author

Sorry, what's the issue?

I'm an idiot. Thanks.

@tlunttil
Copy link
Copy Markdown
Contributor

@tlunttil the goal posts have moved considerably since Reinhold Schaaf's PR was merged in.

Ok, I'd completely missed that somehow. In that case this PR looks good to me now.

@smaret
Copy link
Copy Markdown
Contributor

smaret commented Jul 13, 2016

In that case this PR looks good to me now.

Thanks. @allegroLeiden: could you please merge the master branch into this branch so I can merge your PR?

@allegroLeiden
Copy link
Copy Markdown
Contributor Author

Done.

@smaret smaret merged commit 78ee2c0 into lime-rt:master Jul 13, 2016
@allegroLeiden allegroLeiden deleted the fix_68 branch September 7, 2016 12:49
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

Projects

None yet

Development

Successfully merging this pull request may close these issues.

4 participants