Search code examples
gnuplot

Reading right ascension/declination coordinates in gnuplot

I have a two column file with right ascension/declination coordinates:

18:42:21.8 -23:04:52
20:55:00.8 -17:23:19

I can read the first column specifying data as 'timefmt' but it seems there is no way to do a similar reading for angular data. I could, of course delete :'s and plot ($2+$3/60+$3/3600) but I wonder if there is a more elegant way.

Solution

  • You can define a function which is doing the job for you, which might be a bit more convenient and shorter in the plot command.

    • Convert your hours, minutes, seconds or degrees, minutes, seconds into seconds via strptime() or timecolumn(). In gnuplot console type check help strptime, help timecolumn and help time_specifiers. Use %tH:%tM:%tS, not %H:%M:%S.

    However, you have to be careful how gnuplot interprets negative times:

    if your input time is for example -00:17:56.7 gnuplot will interpret this as +00:17:56.7 which is not what you expect. Apparently, -00 is equal to +00 and hence 17 is interpreted as positive, although you intended it to be negative. A workaround in this special case would be the following:

    Create a function myTimeSign(s) which checks if hours are 0 and if the first character of your time is - and will return -1, and 1 otherwise. 

    myTimeSign(s) = strptime("%tH",s)==0 && s[1:1] eq '-' ? -1 : 1

    Multiply this with your time. This will do here as workaround, but not in general.

    Update: This has been reported as bug (https://sourceforge.net/p/gnuplot/bugs/2245/) and is already fixed in the development version of gnuplot.

    Code:

    ### time / angle conversion
reset session
set size square
set object 1 rect from graph 0,0 to graph 1,1 fc rgb "black"

$Orion <<EOD
05:55:10.29   +07:24:25.3   0.42   Betelgeuse
05:14:32.27   -08:12:05.9   0.18   Rigel
05:25:07.87   +06:20:59.0   1.64   Bellatrix
05:32:00.40   -00:17:56.7   2.20   Mintaka
05:36:12.81   -01:12:06.9   1.69   Alnilam
05:40:45.52   -01:56:33.3   1.88   Alnitak
05:47:45.39   -09:40:10.6   2.07   Saiph
05:35:08.28   +09:56:03.0   3.47   Meissa
EOD

myTimeFmt = "%tH:%tM:%tS"

RA(n) = timecolumn(n,myTimeFmt)
myTimeSign(s) = strptime("%tH",s)==0 && s[1:1] eq '-' ? -1 : 1   # returns -1 if hours are -00
Dec(n) = timecolumn(n,myTimeFmt)*myTimeSign(strcol(n))

set xrange[strptime(myTimeFmt,"06:12"):strptime(myTimeFmt,"05:00")] reverse
set format x "%H^h%M^m" time
set yrange[strptime(myTimeFmt,"-12:00"):strptime(myTimeFmt,"+12:00")]
set format y "%tH°%tM'" time
set tics out

plot $Orion u (RA(1)):(Dec(2)):(-log10($3)+1.5) w p pt 7 ps var lc rgb "yellow" notitle
### end of code

    Result:

    enter image description here