Question about fields on different frames (FK5, ICRS)

[miguelcarcamov]

Hello,

I have just noticed that some datasets might have fields with different frames. One field could be in FK5 and other one could be in ICRS. Is there any taql query of casacore code that allows to retrieve all fields with the same frame using casacore?

Cheers,
Miguel

[kgolap]

When a column in FIELD have varying frame there will be an optional column with an index defining the frame...say PHASE_DIR will have a corresponding PhaseDirRef optional column with an index value ...all rows with the same value will have the same frame...and the frame it associate with will be the frame with the index value in the MEASINFO column keyword of PHASE_DIR

[miguelcarcamov]

When a column in FIELD have varying frame there will be an optional column with an index defining the frame...say PHASE_DIR will have a corresponding PhaseDirRef optional column with an index value ...all rows with the same value will have the same frame...and the frame it associate with will be the frame with the index value in the MEASINFO column keyword of PHASE_DIR

Thank for the response. Exactly, that's what I'm saying. The problem is not how to get the coordinate system, the problem is to have all coodinates on a common coordinate system. I'm asking if that is possible by doing a TaQL query or casacore.

[kgolap]

Well it is as simple as a TaQL to return all the rows with the same value (say 2) in the PhaseDirRef column

[miguelcarcamov]

Well it is as simple as a TaQL to return all the rows with the same value (say 2) in the PhaseDirRef column

Again, you don't understand. The problem is not to know what coordinate system a specific field has (I already know how to do that). I want to know if there is a TaQL query or casacore command that allows me to have a common coordinate system for all my ra and dec coordinates.

[gervandiepen]

Using the meas module in TaQL you can convert the coordinate system of data.
For example:
select meas.j2000(PHASE_DIR) from ~/data/GER.MS::FIELD

Using TaQL's update command you can use that function to change the PHASE_DIR
update ~/data/GER.MS::FIELD set PHASE_DIR=meas.j2000(PHASE_DIR)

When doing this, it is wise to first make a copy of the FIELD table in case the update is not doing what you expect.

Note that the meas module is dynamically loaded by TaQL as the shared libcasa_meas library. I believe that CASA is using static builds, so I'm not sure if it will work properly in CASA.

See also section 4.10.18 (special Measures functions) in TaQL note 199.

[miguelcarcamov]

Hi @gervandiepen , thank you for your response.

Quick question then - since J2000 is just the epoch - let's say I have 2 rows in my FIELD table. One row uses FK5 frame ("J2000" according to measurement set) and the other one uses ICRS. If I do select meas.j2000(PHASE_DIR) from ~/data/GER.MS::FIELD then all RA and DEC for PHASE_DIR are converted automatically to FK5 right?

[gervandiepen]

Hi @miguelcarcamov

I think it handles varying frames correctly, but I'm not fully sure. It's a long time ago this code was written. Underneath it uses TableMeasures which deal with those cases.
You have to try it.

[tammojan]

J2000 in casacore is actually referring to the frame that is close to FK5, so in casacore it refers to both the epoch and the frame simultaneously. I don't think the FK5 frame itself is supported, but to my understanding they are very very close.

[gervandiepen]

@miguelcarcamov Did you manage to change the coordinates?
This morning I realized you also need to change the refcode in the referred column and set it to 0 (meaning J2000).
Alternatively, you can change the column keywords telling the frame is the same for the entire column, but that is a bit harder (although it can be done using TaQL).

[miguelcarcamov]

Hi @gervandiepen, @tammojan I confirm the coordinates change, however, I'm getting a shifted image as result (I am wondering if the J2000 function applies the function to RA and DEC even if they are in J2000). Also I was wondering if casacore has a function to convert from FK5 to ICRS and viceversa, so maybe I can make the conversion after reading the coordinates, and not while reading the coordinates in my taql query.

[gervandiepen]

Hi @miguelcarcamov
Did you adjust the refcode as suggested in my previous message?
You can change the values themselves to any system supported by the Measures system (see note 199). But you also need to change the coordinate system information which requires a separate update command.

[miguelcarcamov]

@gervandiepen I'm not adjusting the refcode since my idea is not to update the measurement set but reading the values correctly, i.e all fields with RA and DEC on the same frame. However, measures system in casacore does not support J2000 -> ICRS conversion, I think I can only convert from ICRS -> J2000. But for me that would be a step back since all people is moving towards working with ICRS. It would be ideal to convert all to ICRS instead of converting all to FK5.

[gervandiepen]

It is possible to convert from J2000 to ICRS. For example in TaQL:

meas.direction('icrs',[30deg,40deg],'j2000')

All systems enum-ed in MDirection.h can be used as to or from. Of course, some also need earth position or epoch.

[miguelcarcamov]

@gervandiepen Then using that function in the TaQL query it should be possible to convert all RA and DEC to icrs. However, I would also need the variable that stores the current field frame. I'm thinking on something like:

select MEAS.DIRECTION('icrs', REFERENCE_DIR, CURRENT_FRAME) as REFERENCE_DIR,MEAS.DIRECTION('icrs', PHASE_DIR, CURRENT_FRAME) as PHASE_DIR,ROWID() A
S ID FROM dataset.ms/FIELD where !FLAG_ROW

However, I don't know how to get the frame keyword from the FIELD table with TaQL - I might need to check if something appears on note 199. I think it has keyword measure, this could also be written as:

select MEAS.DIRECTION('icrs', REFERENCE_DIR) as REFERENCE_DIR,MEAS.DIRECTION('icrs', PHASE_DIR) as PHASE_DIR,ROWID() A
S ID FROM dataset.ms/FIELD where !FLAG_ROW

And I think this latter would work? i.e if the field is already in ICRS it won't do anything. If the field is J2000 will convert it to ICRS, am I right? From the note 199:

If the value is a table column containing measures, the measure frame is obtained from the column and no extra frame argument should be given.

[gervandiepen]

It should be sufficient to do:

select meas.direction('icrs',REFERENCE_DIR) from ~/data/GER.MS::FIELD

The REFERENCE_DIR column knows its frame which is used by the TableMeasure and TaQL.

[miguelcarcamov]

It should be sufficient to do:

select meas.direction('icrs',REFERENCE_DIR) from ~/data/GER.MS::FIELD

The REFERENCE_DIR column knows its frame which is used by the TableMeasure and TaQL.

Awesome! @gervandiepen That's what I thought. Also - could you please give me an example on how to do the same using C++ casacore? I have seen that there is a DirectionEngine class, however I haven't found a proper example. For example, let's say that I have two double variables ra, dec and I know they are framed in J2000. How could I convert them to ICRS?

[gervandiepen]

You can also execute TaQL commands from C++ using function tableCommand in TableParse.h.

You can use the Measures machinery, but that's a bit complex.
I guess you can use DirectionEngine, but I have little time to look at it. You can take a look at DirectionUDF.cc how it is done for TaQL.

[gervandiepen]

I realized this morning that in C++ you can use the TableMeasures to do the conversion. See class ArrayMeasColumn. The MS classes give you such an object for the reference-dir.

[miguelcarcamov]

Hi @gervandiepen thank you very much for your suggestions. Just want to clarify that I have already converted the MS REFERENCE_DIR/PHASE_DIR ra, dec to ICRS using TaQL (yey!). However, I would like to create a more general function that could convert any double ra,dec coordinate (not in a Measurement Set) from J2000 to ICRS. I have seen ArrayMeasColumn but I think the constructor needs a table. But I don't know if that is what I'm looking for. I'm looking for something more like:

double ra = 0.1; //deg - let's say it is in J2000
double dec = 0.2; //deg - let's say it is in J2000
new_ra = convert("icrs", ra, "j2000");
new_dec = convert("icrs", dec, "j2000");

or

double coords[2] = {0.1, 0.2}; // ra,dec in deg in J2000
new_coords = convert("icrs", coords, "j2000");

where convert can be any casacore function that allows me to convert either a double/float or c-array/pointer. I guess what you mean is that I would need to create a table and put those values right?

[gervandiepen]

In that case it's best to use MDirection in the Measures directly.

[miguelcarcamov]

@gervandiepen can you please give me a brief example of how can I use the MDirection class?

[gervandiepen]

See measures/Measures.h and tMDirection.cc

[miguelcarcamov]

@gervandiepen I have done something, however, I'm not sure if it's correct and I don't know how to know on which units the result is, and how to get the results from this. Check this:

My code is:

std::cout << "Original direction is: (" << ra << "," << dec <<") deg" << std::endl;
      casacore::MDirection::Convert conv(casacore::MDirection::Ref(casacore::MDirection::J2000), casacore::MDirection::Ref(casacore::MDirection::ICRS));
      casacore::MDirection coord = casacore::MDirection(casacore::Quantity(ra, "deg"),casacore::Quantity(dec, "deg"), casacore::MDirection::Ref(casacore::MDirection::J2000));
      std::cout << "Converted direction " << conv(coord) << std::endl;
      std::cout << " \n\n" << std::endl;

The result of this is:

Original direction is: (187.706,12.3911) deg
Converted direction Direction: [-0.967885, -0.130965, 0.214584]

Why it returns 3 values? and also on which units are these values, and how can I get them as double variables?

[miguelcarcamov]

@gervandiepen Are these 3 values the matrix values for the conversion? or am I wrong?

[gervandiepen]

The values you see are the direction cosines. The TaQL command

meas.dircos('icrs',[187.706deg,12.3911],'j2000')

gives the same values. You can use MDirection's getAngle function to get the angles.

[miguelcarcamov]

Yep, that made it! Thank you!!!