<div dir="ltr"><br><div class="gmail_extra"><br><br><div class="gmail_quote">On 16 October 2013 22:02, Thomas Robitaille <span dir="ltr"><<a href="mailto:thomas.robitaille@gmail.com" target="_blank">thomas.robitaille@gmail.com</a>></span> wrote:<br>
<blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex">On 16 October 2013 15:58, Maik Riechert <<a href="mailto:maik.riechert@arcor.de">maik.riechert@arcor.de</a>> wrote:<br>

> Hi Tom,<br>
><br>
> That's exactly what I needed. Thanks!<br>
><br>
> It takes around 3s on my machine for a 4200x2800 image. Is that to be<br>
> expected or are there some speedup tricks?<br>
<br>
To be fair, you *are* converting 12 million positions through a WCS<br>
transformation which is a non-trivial transformation so I'm not sure<br>
if there are any speedups if you need all the values. However, the<br>
question is whether you really need the world coordinates for *all*<br>
pixels. What are you using the values for?<br>
<br></blockquote><div><br></div><div><br></div><div style>AST attempts to speed up this sort of thing by using a recursive scheme as follows:</div><div style><br></div><div style><br></div><div style><br></div><div style>
Initial pane = whole array</div><div style>wcs_values = TransformPane( Mapping, pane )</div><div style><br></div><div style><br></div><div style>where the TransformPane function does the following:</div><div style><br></div>
<div style><br></div><div style><br></div><div style>- Do a linear fit between pixel and world coords at a small number of positions spread across the pane. </div><div style>- If the residuals of the fit are small enough:</div>
<div style>   - use the fit to transform all positions in the pane and return them as the function value.</div><div style>- Otherwise, if the depth of recursion is less than some limit:</div><div style>   - divide the pane into four quarters and invoke TransformPane recursively on each quarter pane, gathering the returned wcs values together and returning them to the caller.</div>
<div style>- Otherwise:</div><div style>   - Use the full mapping to transform the complete pane.</div><div style><br></div><div><br></div><div><br></div><div><br></div><div style>This means that for typical mappings that are locally smooth, you are usually using a very simple and fast linear transformation for each point rather than the full might of the complete pixel->wcs mapping. Clearly, you need to set some limit to the depth of recursion - we use the criterion that the number of pixels in the pane must be at least 4 times the number of points needed to do the linear fit. </div>
<div style><br></div><div style>The user specifies the accuracy with which the points are to be found, and this determines the maximum acceptable residuals for the fit. If the residuals are higher, the pane is subdivided.</div>
<div style><br></div><div style>See for instance <a href="http://starlink.jach.hawaii.edu/devdocs/sun211.htx/node430.html">http://starlink.jach.hawaii.edu/devdocs/sun211.htx/node430.html</a></div><div><br></div><div style>
David</div><div><br></div><div> </div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex">
Cheers,<br>
Tom<br>
<br>
><br>
> Maik<br>
><br>
><br>
> On 16/10/13 15:15, Thomas Robitaille wrote:<br>
>><br>
>> Hi Maik,<br>
>><br>
>> Someone asked a very similar question just yesterday on GitHub:<br>
>><br>
>> <a href="https://github.com/astropy/astropy/issues/1587" target="_blank">https://github.com/astropy/astropy/issues/1587</a><br>
>><br>
>> In short, the solution in your case is:<br>
>><br>
>> NAXIS1 = 4256<br>
>> NAXIS2 = 2832<br>
>> x = np.arange(NAXIS1)<br>
>> y = np.arange(NAXIS2)<br>
>> X, Y = np.meshgrid(x, y)<br>
>> ra, dec = wcs.wcs_pix2world(X, Y, 0)<br>
>><br>
>> Any luck?<br>
>><br>
>> Tom<br>
>><br>
>><br>
>><br>
>><br>
>> On 16 October 2013 14:48, Maik Riechert <<a href="mailto:maik.riechert@arcor.de">maik.riechert@arcor.de</a>> wrote:<br>
>>><br>
>>> Hi,<br>
>>><br>
>>> I'm quite new to python and astropy and got stuck at the following which<br>
>>> probably is extremely easy to do.<br>
>>><br>
>>> I created a WCS object out of a FITS header and wanted to get the RA/Dec<br>
>>> coordinates of each image pixel. I tried:<br>
>>><br>
>>> w.wcs_pix2world(np.ndindex(4256, 2832), 0)<br>
>>><br>
>>> where 4256x2832 are the image dimensions. This returns 'ValueError:<br>
>>> object of too small depth for desired array'. I was trying to avoid<br>
>>> allocating an array with all pixel coordinates and thought an iterator<br>
>>> would work too.<br>
>>><br>
>>> What would you recommend for my case (translating every pixel)?<br>
>>><br>
>>> Cheers and thanks,<br>
>>> Maik<br>
>>><br>
>>><br>
>>> _______________________________________________<br>
>>> AstroPy mailing list<br>
>>> <a href="mailto:AstroPy@scipy.org">AstroPy@scipy.org</a><br>
>>> <a href="http://mail.scipy.org/mailman/listinfo/astropy" target="_blank">http://mail.scipy.org/mailman/listinfo/astropy</a><br>
><br>
><br>
_______________________________________________<br>
AstroPy mailing list<br>
<a href="mailto:AstroPy@scipy.org">AstroPy@scipy.org</a><br>
<a href="http://mail.scipy.org/mailman/listinfo/astropy" target="_blank">http://mail.scipy.org/mailman/listinfo/astropy</a><br>
</blockquote></div><br></div></div>