{"action":"created","issue":{"url":"https://api.github.com/repos/GalSim-developers/GalSim/issues/818","repository_url":"https://api.github.com/repos/GalSim-developers/GalSim","labels_url":"https://api.github.com/repos/GalSim-developers/GalSim/issues/818/labels{/name}","comments_url":"https://api.github.com/repos/GalSim-developers/GalSim/issues/818/comments","events_url":"https://api.github.com/repos/GalSim-developers/GalSim/issues/818/events","html_url":"https://github.com/GalSim-developers/GalSim/pull/818","id":184340025,"number":818,"title":"#799 Refactor the middle-layer code to make more of the operations available from Python.","user":{"login":"rmjarvis","id":623887,"avatar_url":"https://avatars.githubusercontent.com/u/623887?v=3","gravatar_id":"","url":"https://api.github.com/users/rmjarvis","html_url":"https://github.com/rmjarvis","followers_url":"https://api.github.com/users/rmjarvis/followers","following_url":"https://api.github.com/users/rmjarvis/following{/other_user}","gists_url":"https://api.github.com/users/rmjarvis/gists{/gist_id}","starred_url":"https://api.github.com/users/rmjarvis/starred{/owner}{/repo}","subscriptions_url":"https://api.github.com/users/rmjarvis/subscriptions","organizations_url":"https://api.github.com/users/rmjarvis/orgs","repos_url":"https://api.github.com/users/rmjarvis/repos","events_url":"https://api.github.com/users/rmjarvis/events{/privacy}","received_events_url":"https://api.github.com/users/rmjarvis/received_events","type":"User","site_admin":false},"labels":[],"state":"open","locked":false,"assignee":null,"assignees":[],"milestone":null,"comments":3,"created_at":"2016-10-20T21:03:04Z","updated_at":"2016-10-25T20:21:45Z","closed_at":null,"pull_request":{"url":"https://api.github.com/repos/GalSim-developers/GalSim/pulls/818","html_url":"https://github.com/GalSim-developers/GalSim/pull/818","diff_url":"https://github.com/GalSim-developers/GalSim/pull/818.diff","patch_url":"https://github.com/GalSim-developers/GalSim/pull/818.patch"},"body":"Issue #799 nominally involved moving the contents of drawShoot, plainDraw, and fourierDraw up to the Python layer so we can more easily incorporate the changes we need for #722.\r\n\r\nI did manage to do that, although I didn't implement all of the new API design that I billboarded at the LSST workshop yet.  I may do more before the hack week, but there's a lot here already, so I think it's worth putting this out as a PR so people can have a chance to look at it before then.\r\n\r\nAs usual, I got a little carried away on the refactoring, so there's a bit more here that I found useful to change in the process of moving these functions to Python.\r\n\r\nNew features:\r\n\r\n* Implemented ImageC as a complex image type.  It can't be written to/read from FITS files, but everything else seems to be working.\r\n* Switched drawKImage to return an ImageC rather than two real images for (re, im).\r\n* Switched InterpolatedKImage to take an ImageC rather than two real images.\r\n* Added the ability to directly take an FFT or inverse FFT of an image without going through the XTable and KTable classes.  This is more efficient, since it doesn't require any extra storage aside from the input and output images.\r\n* Added a `maxSB` method that all GSObjects implement.  When possible it returns the maximum value of the surface brightness at any location in real space.  For all analytic objects and transformations, it's easy to return this exactly.  For compound objects, it can only return an estimate, which is generally an overestimate.  This is fine, since the use case for this is in calculating a good number of photons to shoot, and an overestimate will just use more photons than necessary, which is inefficient, but not inaccurate.\r\n* Added im[x,y] syntax, which seems like we should have had from the beginning, but for whatever reason we never implemented.  im[x,y] = value and value = im[x,y] both work now and are synonyms for `im.setValue(x,y,value)` and either `value = im(x,y)` or `value = im.getValue(x,y)`.\r\n\r\nSome deprecations:\r\n* `wmult` cf. #810 \r\n* `image.at(x,y)` is replaced by `image.getValue(x,y)`.  This is very rarely used, since `im(x,y)` and now `im[x,y]` are more convenient, and the name of the function doesn't really match up with other such methods (e.g `setValue`), so it seems like `getValue` is more intuitive for the rare case when someone needs this method and wouldn't just use the `()` or `[]` syntax.\r\n* The old syntax with `re, im` images for `drawKImage` and `InterpolatedKImage` still work, but give deprecation warnings.\r\n\r\nVarious mostly functionality-preserving middle-layer changes:\r\n* Changed all the `fillXValue` and `fillKValue` functions to directly write to images, rather than `tmv::MatrixView`s.\r\n* Got rid of `fillXGrid` and `fillKGrid`, since they are no longer necessary with the ability to FFT images directly.\r\n* Moved `goodFFTSize` to Image.h and in Python made it available as a classmethod of Image.\r\n* Allow Images to have non-unit step.  This lets the real and imag parts of complex images be views rather than require copying.\r\n* Changed mutating ImageView methods to not be `const`.  This was a confusing notation, which isn't actually required, so I switched it.  Although, it doesn't really make any difference.\r\n* Made the memory allocation for Images align on 16 byte boundaries.  This lets ImageAlloc memory be appropriate for use with FFTW without requiring the use of `FFTW_Array`.  The latter is still around in `FFT.h`, but it and `XTable` and `KTable` are not used for the regular draw commands.  They are still used by SBInterpolatedImage and SBInterpolatedKImage, which I'd like to change (and use the new direct image ffts), but I'll put that off to another issue.\r\n* Changed the way images are wrapped when the size of the FT is smaller than the full size needed according to maxk.  We used to make it Hermitian along the edges, then wrap it, and make it Hermitian again.  It turns out there is a subtle bug in doing this that only shows up at a very low level (and only popped in the tests when I added another Spergel instance with nu>0 to the list of gsobjects in test_flip).  Now we only wrap once which handles both the wrapping and the Hermitianity, and this seems to have fixed the problem.\r\n* Changed the dbg stuff to use a Singleton pattern, which obviates the need for some of the boiler-plate code at the top of many .cpp files.\r\n* Added a bunch of unit tests.  All the new functionality of course, but also added a function called `check_basic` which runs through a bunch of simple tests of flux, maxSB, xValue, kValue, etc. which turned out to be handy in checking my refactoring of the fillXImage and fillKImage functions.\r\n* Fixed a minor bug in how uninitialized Images are handled in a couple places.  e.g. they used to not pickle correctly.\r\n* You can now construct an image from 2D list, like `im = galsim.Image([[1,2,3], [4,5,6], [7,8,9]], dtype=float)`, which saves you the step of converting to a numpy array first.\r\n* You can also now `drawImage` with `nx,ny` or `bounds` and not provide a `scale`, which used to be disallowed.  Now it automatically sets scale to `obj.nyquistScale()`.\r\n"},"comment":{"url":"https://api.github.com/repos/GalSim-developers/GalSim/issues/comments/256164115","html_url":"https://github.com/GalSim-developers/GalSim/pull/818#issuecomment-256164115","issue_url":"https://api.github.com/repos/GalSim-developers/GalSim/issues/818","id":256164115,"user":{"login":"rmjarvis","id":623887,"avatar_url":"https://avatars.githubusercontent.com/u/623887?v=3","gravatar_id":"","url":"https://api.github.com/users/rmjarvis","html_url":"https://github.com/rmjarvis","followers_url":"https://api.github.com/users/rmjarvis/followers","following_url":"https://api.github.com/users/rmjarvis/following{/other_user}","gists_url":"https://api.github.com/users/rmjarvis/gists{/gist_id}","starred_url":"https://api.github.com/users/rmjarvis/starred{/owner}{/repo}","subscriptions_url":"https://api.github.com/users/rmjarvis/subscriptions","organizations_url":"https://api.github.com/users/rmjarvis/orgs","repos_url":"https://api.github.com/users/rmjarvis/repos","events_url":"https://api.github.com/users/rmjarvis/events{/privacy}","received_events_url":"https://api.github.com/users/rmjarvis/received_events","type":"User","site_admin":false},"created_at":"2016-10-25T20:21:45Z","updated_at":"2016-10-25T20:21:45Z","body":"> How will one access the complex array at the C level?\r\n\r\nIn the C++ layer, it's an array of `std::complex<double>`.  That's just two doubles sequentially in memory: real then imag.  That's the same way numpy stores things for an array of `np.complex128`.  It's also the way FFTW expects things for their complex type `fftw_complex`.\r\n\r\n[Here](https://github.com/GalSim-developers/GalSim/blob/%23799/pysrc/NumpyHelper.h#L202) is the line that casts the C++ memory into a complex numpy array.  T = `std::complex<double>`, and `getCode()` will return `NPY_COMPLEX128`.  The [line](https://github.com/GalSim-developers/GalSim/blob/%23799/pysrc/NumpyHelper.h#L122) to go back the other way from numpy data to `std::complex<double>*` is even simpler, since it's just a `reinterpret_cast`."}}
{"actor":{"display_login":"rmjarvis"}}