Map to FITS file question

206 views
Skip to first unread message

Trae Winter

unread,
Apr 26, 2012, 4:09:32 PM4/26/12
to su...@googlegroups.com
I am manipulating AIA FITS files in Python using the SunPy map object and associated methods.  I need to write out the map object in a format that IDL can read.  I would prefer to save it as a FITS file since that is a "universal" format.

Is there a way to easily write out a SunPy map object to a FITS file?

I've started going through the PyFITS documentation to find an answer, but I thought that someone in the group might know off the top of your heads.

Thanks,
-Trae

Keith Hughitt

unread,
Apr 27, 2012, 10:16:30 AM4/27/12
to su...@googlegroups.com
Hey Trae,

At the moment there isn't any really simple way to save a SunPy Map as a FITS file, although that is something that we should have.

I'll see if I can add a save() method real quick.

Incidentally - do you have any thoughts on how we should handle the FITS headers in SunPy maps? For now, we decided to only store the original header. This way we don't have to update things in multiple places when the map is change (e.g. map.fits_header['cdelt1'] and map.scale['x']). It should not be too difficult to write back all of the necessary parameters (scale, centering, etc) into a new header when saving the map to a file, but do you think it would be very useful to have a copy of the header that is always up to date?

Keith


--
You received this message because you are subscribed to the Google Groups "SunPy" group.
To post to this group, send email to su...@googlegroups.com.
To unsubscribe from this group, send email to sunpy+un...@googlegroups.com.
For more options, visit this group at http://groups.google.com/group/sunpy?hl=en.

Joe Hourcle

unread,
Apr 27, 2012, 10:20:08 AM4/27/12
to su...@googlegroups.com

On Apr 27, 2012, at 10:16 AM, Keith Hughitt wrote:

> Hey Trae,
>
> At the moment there isn't any really simple way to save a SunPy Map as a
> FITS file, although that is something that we should have.
>
> I'll see if I can add a save() method real quick.
>
> Incidentally - do you have any thoughts on how we should handle the FITS
> headers in SunPy maps? For now, we decided to only store the original
> header. This way we don't have to update things in multiple places when the
> map is change (e.g. map.fits_header['cdelt1'] and map.scale['x']). It
> should not be too difficult to write back all of the necessary parameters
> (scale, centering, etc) into a new header when saving the map to a file,
> but do you think it would be very useful to have a copy of the header that
> is always up to date?

I'd suggest doing lazy updates -- give an API to access the header, and
if someone uses it, make sure everything's up-to-date before giving it back
to them / saving the file / whatever.

-Joe


Keith Hughitt

unread,
Apr 27, 2012, 10:25:44 AM4/27/12
to su...@googlegroups.com
That is probably not a bad idea. Originally I thought of there being three different approaches:

(1) Only store the original header
(2) Only store an up-to-date header
(3) Store both the original header and an up-to-date one

Having something like map.get_header() to generate an updated header would probably be useful though and could be used for a .save() function.

This way the user could do something like:

x = sunpy.make_map('file.fits')

...hack hack hack...

x.fits_header           // => original header
x.get_header()          // => updated header
x.save('new_file.fits') // => new fits file with updated header/data


The only downside to this approach is that the way you access the original/updated header is different (property vs. method), but it does save us from having to update properties in multiple places whenever
something is changed.

Keith

Joe Hourcle

unread,
Apr 27, 2012, 10:43:18 AM4/27/12
to su...@googlegroups.com
The other option is:

x.fits_header // original values
x.update_fits_header()
x.fits_header // current values

But that's got its own problems. (and I see no problem with not
updating 'fits_header', as when it's in python, it's not in FITS)

I'd suggest making it a method for *all* access, and deprecate
the direct property access, so that you have the freedom to
change the underlying storage in the future.

-Joe

Keith Hughitt

unread,
Apr 27, 2012, 12:14:34 PM4/27/12
to su...@googlegroups.com
Okay, I submitted a pull request which adds a "save" method to the Map objects. If you have the chance, take a look and let me know what you think and if you have any suggestions.

Keith


-Joe

Trae Winter

unread,
Apr 27, 2012, 12:19:25 PM4/27/12
to su...@googlegroups.com
I'm not sure if I see the value of maintaining the original FITS file header information. Presumably, you would keep the original FITS file that was used to make the map object and could always read it back in if you needed the unadulterated header.  I understand not wanting to update the same information in multiple places, however, I think that writing out an updated file in a format that is agnostic to programming language is an important feature. Sorry my Python skills aren't quite up to the task yet, but I am learning quickly. (I still hate object oriented programming, but I'm sure I'll come around.  Eventually.)

Headers only take up a small amount of storage space, so if the consensus is to keep them both, that's fine.  The only thing I worry about is confusion.  People used to map objects are used to the structure tags automatically updating.  In the save property I would automatically have the saved headers updated.  Half the skill in programming is human engineering and saving users from ID10T errors.   

-Trae

On Fri, Apr 27, 2012 at 10:43 AM, Joe Hourcle <one...@annoying.org> wrote:
--
You received this message because you are subscribed to the Google Groups "SunPy" group.
To post to this group, send email to su...@googlegroups.com.
To unsubscribe from this group, send email to sunpy+un...@googlegroups.com.
For more options, visit this group at http://groups.google.com/group/sunpy?hl=en.




--

"Be yourself, everyone else is already taken." 

- Oscar Wilde


DVD PS

unread,
Apr 27, 2012, 12:34:29 PM4/27/12
to su...@googlegroups.com
 
I'm not sure if I see the value of maintaining the original FITS file header information. Presumably, you would keep the original FITS file that was used to make the map object and could always read it back in if you needed the unadulterated header.  I understand not wanting to update the same information in multiple places, however, I think that writing out an updated file in a format that is agnostic to programming language is an important feature. Sorry my Python skills aren't quite up to the task yet, but I am learning quickly. (I still hate object oriented programming, but I'm sure I'll come around.  Eventually.)

Headers only take up a small amount of storage space, so if the consensus is to keep them both, that's fine.  The only thing I worry about is confusion.  People used to map objects are used to the structure tags automatically updating.  In the save property I would automatically have the saved headers updated.  Half the skill in programming is human engineering and saving users from ID10T errors.   

Wow Keith! that's was quick!
I'm just would like to suggest to have in the COMMENTS field information of when, and what's been used to saved the file, so the people could differentiate between the original and the modified file (I'm sure some people uses same file names!)

my 2 cents ;-)

David

Ken Dere

unread,
Apr 28, 2012, 11:14:55 AM4/28/12
to su...@googlegroups.com
In my opinion, Hierarchical Data Format (5) is a better way to go.  It is well supported in Python with pytables and has excellent support in IDL.  FITS is not officially supported by IDL (last I looked).  The original FITS file format was to serve a single purpose of transporting images.  The binary FITS file is a Frankenstein creation.

Ken Dere

Keith Hughitt

unread,
Apr 29, 2012, 8:11:28 AM4/29/12
to su...@googlegroups.com
Hi Ken,

It looks like there is an even lighter-weight Python binding which is almost exactly the same as the original C library (http://code.google.com/p/h5py/) which seems to play well with NumPy.

It should be pretty easy to add support in sunpy.io and sunpy.map.save for reading/writing to HDF5. Any chance you would like to help add support? :)

Best,
Keith


--
You received this message because you are subscribed to the Google Groups "SunPy" group.
To view this discussion on the web visit https://groups.google.com/d/msg/sunpy/-/xcXka0ehqEEJ.

Keith Hughitt

unread,
Apr 29, 2012, 8:12:56 AM4/29/12
to su...@googlegroups.com
Hey David,

I thought about that as well. Do you think we should just append some string (e.g. saved using SunPy xxx on yyy) to the end of the current comment string?

Keith

--

DVD PS

unread,
Apr 30, 2012, 5:44:23 AM4/30/12
to su...@googlegroups.com
On 29 April 2012 13:12, Keith Hughitt <keith....@gmail.com> wrote:
Hey David,

I thought about that as well. Do you think we should just append some string (e.g. saved using SunPy xxx on yyy) to the end of the current comment string?

I think IDL just adds this to the SIMPLE keyword:
SIMPLE  =                    T /Written by IDL:  Mon Apr 30 08:06:52 2012      

I don't really understand why that's written there... maybe to show it as the first thing you see when looking at the header (??). In our case we could provide the original filename in the Comments, plus what you just said, Sunpy version blah, user blah date blah...

David

PS. I don't know if there are standard keywords for that in the FITS description. I will look it up and let you know what I find.

Ken Dere

unread,
Apr 30, 2012, 9:03:23 AM4/30/12
to su...@googlegroups.com
Hi Keith,

I really just do not have the time to work on that.

best

Ken
To unsubscribe from this group, send email to sunpy+unsubscribe@googlegroups.com.

Joe Hourcle

unread,
Apr 30, 2012, 10:35:40 AM4/30/12
to su...@googlegroups.com

On Apr 29, 2012, at 8:12 AM, Keith Hughitt wrote:

> Hey David,
>
> I thought about that as well. Do you think we should just append some
> string (e.g. saved using SunPy xxx on yyy) to the end of the current
> comment string?

So, as I actually analyzed a bunch of FITS files ... IDL is the
*only* one that writes it at the top of the file. (and of course,
as IDL doesn't have direct support, it was obviously written by
a specific package in IDL, but it never tells you)

And that often means that it masks the fact that it's a FITS file,
as it never clearly says 'this is a FITS file' anywhere in the file.

Some programs use 'ORIGIN' to store the software information,
but that's actually a mistake ... ORIGIN is supposed to be for
the institution that created the file, not the software that
created the file.

(I won't say where these came from .. but the second one is NOT
from NOAO)
ORIGIN = 'wrfits ' /
ORIGIN = 'NOAO-IRAF FITS Image Kernel July 2003' / FITS file originator


I'd say that the best place to put it would be a HISTORY card
(which is behaves like COMMENT for the most part, but it's
specially labeled)

eg, if we take a look at a STEREO/SECCHI/EUVI file:


HISTORY Id: make_scc_hdr.pro,v 1.169 2011/10/12 19:06:56 nathan Exp
HISTORY Id: getsccsecpix.pro,v 1.15 2008/05/02 17:09:59 nathan Exp
HISTORY Id: getscccrpix.pro,v 1.10 2011/09/19 22:12:23 nathan Exp
HISTORY Id: getsccpointing.pro,v 1.11 2011/09/21 21:23:43 nathan Exp
...
HISTORY Id: secchi_reduce.pro,v 1.386 2011/09/15 21:57:35 nathan Exp
HISTORY Processed on iapetus: x86 linux IDL7.0
HISTORY Id: scc_icerdiv2.pro,v 1.19 2011/09/15 21:57:35 nathan Exp
HISTORY Id: secchi_rectify.pro,v 1.22 2009/12/14 18:49:22 nathan Exp
HISTORY Id: euvi_point.pro,v 1.14 2011/02/02 12:54:32 mcnutt Exp

although, they have nothing in there about what program wrote it out,
other than the horrible

SIMPLE = T / Written by IDL: Fri Nov 4 12:12:24 2011

Others have stuff in HISTORY such as:

HISTORY Created by tiftofits version 1 on 2006-02-08 at 22:56:40

HISTORY created by MDI_image_reconstructor Wed Jan 6 13:26:02 2010

HISTORY V15 18 Jun 1997 MAKE_FITS_HDR

HISTORY 'KPNO-IRAF' /

(yes, that last one's strange ... technically, the whole thing is a comment)

GONG and SDO/HMI use a numbered field to store the processing details, but I'm
really not a fan ... GONG at least has another field w/ the count of how
many fields to look for, but I don't really see any benefit of that over
just using HISTORY.

-Joe


ps. HDF5 might be better for some types of data (maybe even all
types), but getting solar physicists to give up both IDL *and* FITS?
That's just talking crazy. If you're trying to support other
disciplines who also need solar images (magnetospheric folks, radiation
belts, ITM), I could see trying to support CDF, HDF and NetCDF, but I
assume you'd want to know how many people fit those groups and are using
SunPy to determine how to devote to it. (in part, as you'd need people
who know what a real HDF / NetCDF / CDF file is supposed to look like).





Keith Hughitt

unread,
Apr 30, 2012, 10:53:53 AM4/30/12
to su...@googlegroups.com
Hey Joe,

What do you think an "ideal" history comment from SunPy would look like, and what information would it include?

Keith

DVD PS

unread,
Apr 30, 2012, 9:49:56 PM4/30/12
to su...@googlegroups.com
Hey Joe!

Your description is quite amazing and very precise.  Nevertheless I've just finished to go through the FITS documentation (http://fits.gsfc.nasa.gov/fits_standard.html) to see if there was something else to add.  And nope, not at all.

I completely agree with you that we should use the HISTORY keyword to add what's new on the image... Now, I don't know how easy is to do so in python... as probably there's not an easy way to record every single operation that's done on the array... (or maybe it is).
So, an ideal-world example would be:

HISTORY Data processed by "name" using SunPy v1.2 yyyy-mm-ddThh:mm:ss
HISTORY Image aligned map.align('map-1')
HISTORY Image clipped map.clip(0,100)
HISTORY Image rotated map.drot(3h,'method a')
HISTORY Image smoothed map.smooth(5)
etc etc...

So, this would allow to anyone to see what's done in the image.  Now how useful this will be?  I think in a way will make the people to trust more in the processed data, and share that one rather than the original, the software and then to re-process everything...

David

Keith Hughitt

unread,
May 1, 2012, 10:14:36 AM5/1/12
to su...@googlegroups.com
I was thinking about that. What we could do is create a "history" module which tracks operations performed on a file. Basically for each method called, the method would call something like history.record('....'). Then when it comes time to write the file out, save() could call "history.get_history()" to get a string with all of the changes made.

It probably wouldn't be difficult to record an exact replica of the steps taken (e.g. for "map1 - map2" you don't want to store the entire map2 data in the history of map1), but we could still include the basic operation.

If we do include this kind of providence information, should this be required or optional when saving a file, and if it is optional, should it be on by default?

For the user name, we could include fields in the newly-added .sunpyrc file for name and email address. If they are not set we could either prompt the user, leave the user info out, or allow the user to specify that info via keyword args.

Best,
Keith

Reply all
Reply to author
Forward
0 new messages