BenningtonCS / Telescope-2014

4 stars 0 forks source link

Figure out what the heck everything is in srtn #118

Closed theDarkLard closed 9 years ago

theDarkLard commented 9 years ago

what's going on in there?? Displays, small words with numbers next to them, record files oh my

edaniszewski commented 9 years ago

If ya need any help with analyzing the code to figure out what it does, I can try to help.

I added file header comments to some files describing what they did, but for the more nebulous, poorly documented things.. I'm around.

theDarkLard commented 9 years ago

So the spectra that are covered by the old srtn documentation, i.e. documentation of the java version, only cover two different spectra: discrete spectra, and accumulated spectra. The accumulated spectra is definitely present in our srtn; it is the bottom spectra and denoted by "av spectrum". This is what the old documentation has to say about the accumulated spectrum:

"The [bottom] spectral-line plot shows the accumulated spectra since the selected observation began. Listed at the top of this window is the title "av. Spectrum" as well as the total "integration" time. The bottom script shows the same max-min difference and increasing frequency direction." We don't have a bottom script and our srtn only displays one integration time, which i am guessing is the total...?

I believe that the top two spectra are two different discrete spectra...but note the verb "believe." Here's what the old documentation says:

"The [top..?] is the plot of each individual spectrum as it finishes the user input span in MHz (see "freq" in section 1.3.1). The top of the display lists the input center frequency and the frequency step. The bottom line lists the difference of the maximum and minimum values measured during the frequency scan. There is also an arrow indicator showing the direction of increasing frequency."

With this description in mind, let us revisit all three spectra. The top spectra says 2.4 MHz at its top, which I believe is the frequency span. So for whatever the user sets as the center frequency in freq, the telescope will scan 1.2 MHz above and below this. ON the bottom we have Tsys (which is obvious, and unchanging as long as "cal" is not run) and smax. Still not quite sure what smax is. Now on the the middle spectra. At the top, there reads 0.2-2.2 Mhz, which i think is another frequency span (again, speculation). On the bottom, we have "fs 67.15K pwr 1566.2" "these values change every second). I think fs refers to that difference b/t max and min values during the frequency scan, since there is a K beside it. And pwr of course means power, but I am still unsure of what units those data are in. Will be digging around in the software to try and find out. The bottom spectrum is the sam as the middle, except with "av. spectrum" at the top, as I said before.

theDarkLard commented 9 years ago

@edaniszewski yeah i will probably need some help. Will let you know. :)

theDarkLard commented 9 years ago

This is what the MIT SRTN manual has to say about the output files:

The output data file (yydddhh.rad) is an ASCII text file. Data reduction on the raw output can be done with a spreadsheet program like MS Excel with some effort. The development of spreadsheet MACROS to reduce the data are desirable for large files and long integrations. Future releases or updates of the SRT software could also output certain observations in formats for other data reduction methods. Comment lines: Start with an asterisk and list the STATION LAT and LONG(E/W) on the first line.

So i'm sure we can write a script that pinpoints the relevant data and extracts it. There is also a script to do this in matlab already made, via http://www.astro.gla.ac.uk/observatory/srt under the heading "Analyzing Data".

And the document that I've been pulling all this info from can be found here: http://www.haystack.mit.edu/edu/undergrad/srt/SRT%20Software/SRTManual.pdf

bgwalter commented 9 years ago

I'll also put myself out there as someone who can help figure stuff out

theDarkLard commented 9 years ago

yay friends!!

hcrowl commented 9 years ago

Who needs other social networks when there's github!

theDarkLard commented 9 years ago

Ok so I was not able to come up with much with my limited computer skills. Things we don't know for sure what they are yet are: fs, smax, pwr, and what their units are. I found these variables in the files main.c, plot.c, velspec.c, vspectra.c and the rest of the vspectra series. Maybe those of you with greater facility in computers can put the pieces together and see how these variables are working throughout all the files? Other useful variables might be Tant, which refers to antenna temperature, and Tsys, which refers to system temperature. I should probably learn c :/

edaniszewski commented 9 years ago

Github is the best social media!! :octocat: :bangbang:

Things we don't know for sure what they are yet are: fs, smax, pwr, and what their units are.

@theDarkLard I'll probably look a little bit into the offending variables to try and figure out what their deal is later today after some sushi and Sopranos. (gotta have priorities)

There is also a script to do this in matlab already made, via http://www.astro.gla.ac.uk/observatory/srt under the heading "Analyzing Data".

Good find. Looks like they have some useful resources there, but it also looks like they may be using an older version of the software, so keep that in mind. Using the matlab script would be good, if you have matlab. I hear its pretty expensive to license (though maybe cheaper if the school licenses it?). Rolling out our own allows us to have more control and flexibility with it though. There could be options to save it out as CSV, XML, JSON, etc. so different plotting/analysis modules could be built on top of it. It would also allow it to be exposed via a REST API (assuming things get complex later on).

I should probably learn c :/

I'm an advocate, but I'm not sure that this software need be written in C. Yes, C is fast, but Java can be almost as fast in some cases(http://benchmarksgame.alioth.debian.org/u64q/java.html), and its easier to read (though it has some problems of its own). Python/Cython could even be used, with the computationally intensive things compiled in C, and called by the Python code. I've kind of been wanting to do something like this for a while, but I'm not really sure how to start/if I would remain sane. Regardless, UI work in almost any language seems easier than it is in pure C.

edaniszewski commented 9 years ago

Things we don't know for sure what they are yet are: fs, smax, pwr, and what their units are.

fs Not too sure what this is. It looks like the only place this is mentioned is when it is printed to the UI. Makes me think that is an abbreviation for something well known in radioastronomy. @hcrowl any ideas? It could be describing a temperature since it seems to have units K (Kelvin?) associated with it.

Relevant code snippet:

if (!d1.bsw || iav == 1) {
    if (d1.caldone)
        sprintf(txt, "fs %5.2fK Tant %4.1fK", ddt, d1.tant);
    else{
        printf("fs %5.2fK pwr %4.1f\n", ddt, pwr);
        snprintf(txt, 80, "fs %5.2fK pwr %4.1f", ddt, pwr);
    }
} else
    sprintf(txt, "fs %5.2fK bswpwr %5.2f", ddt, d1.bswpwr);

smax I have no idea about this one. I've looked at all occurrences in the code and was not really able to determine what it was, since everything else it was associated with was just as, if not more, nebulously defined.

Its clear that it is the max value for something represented by 's'. What that 's' is, I am not sure.. From the code snippet below, where it prints out a message, I am making a 100% guess based on its proximity to words starting with 's' that it could be the max sigma value? Take what I say with a barrel of salt though

if (max > spec[maxi / m] + d1.rfisigma * noise && d1.printout) // rfisigma sigma
    printf("check for RFI at %8.4f MHz max %5.0e av %5.0e smax %5.0f %3.0f sigma\n", maxi * d1.bw / blsiz2 
          + d1.lofreq, max, spec[maxi / m], smax, (max - spec[maxi / m]) / noise);
}

pwr This one is interesting. I've seen it used nebulously a lot in the code. Perhaps it is used in different ways. One revealing bit:

d1.tant = pwr - d1.tsys;

that is to say, the temperature of the antenna is equal to the power minus the temperature of the system. Could the units of power also be Kelvin? Or is this just a poor variable naming choice?

All other occurrences of pwr I've seen are pretty nebulous. Perhaps its worth contacting Alan Rogers to get a better understanding of the variable names?

edaniszewski commented 9 years ago

From srt.hlp:

cal

Calibration - various modes are available enter desired mode in srt.cat file using CALMODE keyword. Suggest CALMODE 2 or 3 The basic method is to get the system temperature from the "Y-factor"

  Y-factor = power_on_hot_load / power_on_cold_load
           = (tsys + hot_temperature) / (tsys + cold_temperature)

  solving   tsys = (hot_temperature - y*cold_temperature) / (y-1)

power_on_hot_load = with absorber over the feed or pointing at the trees power_on_cold_load = looking at cold sky hot_temperature = ambient ~ 290 K cold_temperature = 3 K from CMB

Since the LNA noise temperature is very stable is in not really necessary to get the system temperature every time you point the telescope. What is needed is to get the power looking at an absorber or pointing at the trees fairly frequently and the the gain of the LNA + amplifiers + dongle change with time. CALMODE 3 is appropriate if you can point at the trees. In this case just click on cal and the software will divide by a reference spectrum and convert the power to kelvin assuming the power while still pointed at the trees corresponds to Tsys + Tcal. CALMODE 2 does the same by waits for the operator to hit enter when the absorber has been placed over the feed.

also:

mode

The SRT operates in a total power mode. The bandpass is normalized by clicking the cal botton. NPOINT does 5x5 point scan to find peak power and pointing offsets. BEAMSW can be used for sources which have suitable reference spectra on either side of the source. The beam throw equals the BEAMWIDTH set in srt.cat. NBSW sets the number of integration periods in each beam position.

possibly this helps describe what pwr may mean?

theDarkLard commented 9 years ago

Yeah I think that's part of it ^^ I poured through the help file yesterday, all this stuff is up on the wiki. That clip of code you threw up there is interesting though...I saw it but didn't put the pieces together. This is the best evidence yet for pwr being in Kelvin, and it tells us where the numbers for pwr are coming from...pwr is starting to make a little more sense, thanks a ton Erick! Also here's a little Kelvin to Watt conversion which can then be converted to Jansky's:

Noise temperature in Kelvin. Noise temperature in Kelvin = Watts/1.38E-23/1Hz. 1.8E+9 K = 25.1E-15 W/1.38E-23/1 Hz. Noise temperature is in per bandwidth of 1Hz.

I got this from http://www.stargazing.net/david/radio/UnitsConversion.html

hcrowl commented 9 years ago

A question for @theDarkLard:

What does it mean to measure a "power" (ie a flux) in Kelvin? This is a fundamental unit mismatch that shouldn't work. So, what assumption is being made to allow this to happen.

(This is something that I think I know the answer to, but would have to do the reading to make sure. Why don't you do the reading and write a little something up?)

theDarkLard commented 9 years ago

OK so I pulled this info from this webpage: https://casper.berkeley.edu/astrobaki/index.php/Single_Dish_Basics#Some_Important_Equations

"Recall that we relate the power received by the antenna at a frequency ν, Pν, to the antenna temperature, TA, and the observed intensity, Iν, as follows: ed55b78ad411cb36e70217c3ea29dc35 Therefore, the conversion between K and Jy is just fdd3d018fb00dcfc155b9e36068430fa "

We know our frequency and are getting constant antenna temp measurements. Ae is effective area of the antenna, which can be calculated from "Ae = ηaAp (ηa is the aperture efficiency and Ap is the projected area of the telescope)," which can both be tested for our telescope. I don't think we need to do this, just figured it would be good background for where this conversion is coming from.

If this is not exactly what @hcrowl had in mind, I think the other important formula is something called the Johnson-Nyquist formula, P = kTdv where P is power going into the receiver (dongle), k is Boltzmann's constant, and T is antenna temperature. Since k is a constant and dv is just a differential perhaps the srtn software is compartmentalizing those terms and going straight for the P=T, which in our case is P = Tsys + Tant. That's just the logic I am using to try to justify what's in the code.

theDarkLard commented 9 years ago

And to expand a bit, I think the v from dv is actually a greek "noo"...? which is standing in for frequency, so it's like power for that specific frequency..? The flux density equation will help this make more sense, which is what we want to measure anyway. Flux density is how much energy is passing through a surface at a given time. "The Flux Density is equal to the Planck function = monochromatic intensity - I(v) integrated over solid angle of telescope): screen shot 2015-06-23 at 11 15 05 am i.e. the Flux density is just the brightness temperature integrated over the source."

Got this info from https://www.astron.nl/~mag/dokuwiki/lib/exe/fetch.php?media=radio_astronomy_lec_3_ma_garrett.pdf

theDarkLard commented 9 years ago

And notice the units, W/m^2/Hz. Power per unit area per frequency, so it's assumed that you must know what frequency your antenna is tuned to. (i think)

edaniszewski commented 9 years ago

This is the best evidence yet for pwr being in Kelvin, and it tells us where the numbers for pwr are coming from...

I agree that the units for pwr are in Kelvin at some point during the processing, but I wouldn't say with absolute certainty that the pwr with units Kelvin is what is being displayed. It could just be used for some internal calculations. Its kind of tough to follow where the actual conversion occurs (and there is a conversion, as explicitly stated in the help file), and from what units it is converted from.

Having the equations help though, since they would presumably be somewhere in the code.

theDarkLard commented 9 years ago

I think we're good on this now.