Mad-Science
  • Latest Posts
  • Astronomy
  • Geology
  • Natural History
  • Downloads

Radio Astronomy

  • Magnetometers
  • VLF Receivers
  • Muons
  • H-Line Radio Telescope
  • Live Instruments

Muon Counting - a brief update

Parent Category: Radio Astronomy
Published: 20 November 2024

I've just spent a couple of evenings reviewing the code that displays the Muon Counter detail. One of the problems of 'borrowing' other peoples code from the internet is that you may not fully understand it when you need to edit the code. My unhappiness stemmed from the fact that the existing graphs were not really showing me anything, I think this is because the sample window (was 1 minute) was too long and I was losing valuable data in the general noise.

The graphical display code has been re-written to read blocks of data from the CSV file in a 10 second window. I can change this fairly easily in the future as I have actually commented my code.

The revised graphs are now very dense (8640 display points to the graph) and so I am not sure how this will work with a graph width of 1600 pixels.

I may need to rethink.

For information (if you are suffering from insomnia), the process is a broadly as follows:

1. I read a daily block of data from the CW Muon Counter Raspberry Pi.

It looks like this truncated file from 2024-11-17: 

Y,M,D,H,M,S,Event, Ardn_time[ms], ADC[0-1023], SiPM[mV], Deadtime[ms], Temp[C], Name
2024,11,17,00,00,11.869169,4959769, 1032619463, 453, 113.81, 845053568, 23.15,WOBO_MUON_1
2024,11,17,00,00,32.307766,4959770, 1032639923, 496, 144.49, 845057408, 23.15,WOBO_MUON_1
2024,11,17,00,00,33.598346,4959771, 1032641193, 377, 76.28, 845057600, 23.15,WOBO_MUON_1
2024,11,17,00,00,34.803068,4959772, 1032642419, 389, 81.58, 845057984, 23.15,WOBO_MUON_1
2024,11,17,00,00,44.020615,4959773, 1032651632, 415, 95.00, 845059712, 23.15,WOBO_MUON_1
2024,11,17,00,00,46.323831,4959774, 1032653933, 369, 73.42, 845060096, 23.15,WOBO_MUON_1
2024,11,17,00,00,53.045375,4959775, 1032660652, 425, 99.48, 845061440, 23.15,WOBO_MUON_1
2024,11,17,00,01,03.371836,4959776, 1032670971, 728, 431.28, 845063360, 23.15,WOBO_MUON_1
2024,11,17,00,01,15.304950,4959777, 1032682898, 472, 125.87, 845065664, 23.15,WOBO_MUON_1
2024,11,17,00,01,19.087461,4959778, 1032686678, 704, 326.52, 845066432, 23.15,WOBO_MUON_1
2024,11,17,00,01,24.596963,4959779, 1032692184, 450, 114.31, 845067392, 23.15,WOBO_MUON_1
2024,11,17,00,01,29.826490,4959780, 1032697220, 352, 66.89, 845068352, 23.15,WOBO_MUON_1

This is changed from the standard muon data file in that all the commas (its a comma separated file) are written at source by modifying the 'C' code on the CW counters. The commas delimit the tabular columns so that the individual date and time data elements (Y,M,D,H,M,S) can be easily separated for post processing.

 

To produce the graphs, the only data I need is the time (H,M,S) for each detection.

I then create a separate array and populate the first column of the table with the number of seconds in a day divided by the sample window, currently 10 seconds. 

The array has bounds of 8640 rows x 2 columns.

I simply then parse through the imported csv file, calculate the time in seconds for each event, divide that by the sample window value and populate the array. In the sample above, the three highlighted entries would appear as a count of 3 for the sample window starting at 00:00:30 and ending at 00:00:39.

The array is then passed to a Python module (MatplotLib) to produce the graph.

One problem that I haven't as yet been able to solve is how to set the X axis values from a weird seconds /10 value to simple hours. 

 

Postscript:

After reviewing the following days results, I could see that many of the bars were not being plotted. 

The program was then modified to plot 8 x 3 hour windows which is much more of a handful but all the data is visible. I have also annotated the maximum daily count below the 'x' axis to help identify any interesting spikes.

 

Muon Data - Daily Report

Parent Category: Radio Astronomy
Published: 08 July 2024

Muon Data Daily Report

24 hour view

Detailed view

00:00 to 03:00 UTC

03:00-06:00 UTC

06:00-09:00 UTC

09:00-12:00 UTC

12:00-15:00 UTC

15:00-18:00 UTC

18:00-21:00 UTC

21:00- 24:00 UTC

 

 

Link to last 7 days of Muon Data: Muon Counts - 7 day report

Muon Data - 7 day summary graphs

Parent Category: Radio Astronomy
Published: 27 May 2024

The last 7 days of Counting Muons. Not a lot happens here. My average count is approximately 0.242 Muons per second passing through the detectors including other spurious (but rare) events. I'll pick out any 'WOW!' events and post separately on those.

 

Caption

 

 

 

 

 

 

 

Muon Detector - Happy birthday

Parent Category: Radio Astronomy
Published: 10 March 2024

Today, 2024-10-03, is my Muon Detectors first birthday. It operates in co-incidence mode which mean that a count is only registered when a particle passes through both detectors at a relativistic velocity. Hopefully a muon.

Although I have been monitoring the daily count for 12 months, nothing particularly untoward has occurred. To be honest, I wasn't expecting very much, sometimes, patience has to be a virtue. The only gap was caused by a power outage which took a while to spot.

So, I decided to look at the daily counts over the past 12 months and was rather surprised to see a change in the pattern of the count which occurred about 6 months ago. Nothing has changed, the detectors are still in the same place, plugged into the same Raspberry Pi. They have not moved at all.

I am really puzzled as to the cause of this apparent cyclic variation.  I have double checked the data and it appears correct. The average count is 0.18 hits/second.

Investigation ongoing!

UKRAA Cosmic Ray Muon Detector

Parent Category: Radio Astronomy
Published: 11 March 2023

Background

I have recently completed a pair of UKRAA Cosmic Ray Muon Detectors. These are identical to the MIT CosmicWatch detectors and detailed here: MIT

As usual for me, it wasn't quite straightforward. I chose to build from kits and the first unit completed successfully. It then failed after an hour when the SiPM (silicon photomultiplier) failed. This was the most expensive item in the whole kit and I was alarmed to see that there was a 14 month wait time for a replacement (from Farnell). However, this was drastically reduced and I received my replacement after only a 3 month wait. The second unit completed OK without any problems.

 

Two units are required, one operates in 'Master' mode and detects all Muons plus any background radiation. A count is generated for each interaction (hit). Each hit triggers a signal that is passed to the Slave unit. The second 'Slave' Units also detects all Muons and Background Radiation. However, a count is only triggered when there is direction correlation with the trigger signal from the Master at (almost - within 30μS) the same time that the Slave unit detects a hit. There is a small chance that both units could simultaneously detect a hit from a different source - e.g. from background radiation, but this will not make a significant impact on the overall count.  

What are Muons?

For a full explanation read this Wikipedia Entry.

Muons are fundamental charged particles with a mass significantly higher than an electron. They are unstable and have a half-life of 2.2μS. They are formed when there is a high energy interactions with normal matter. In astronomical terms, they are formed when high energy Cosmic Rays collide with atoms (typically Oxygen or Nitrogen) in the upper atmosphere. Most Cosmic Rays are detected from our own Galaxy and a burst could be indicative of a Supernova explosion or activity from our own Sun. These decay down to Muons which are relativistic, travelling at about 99.8% of the speed of light. The Muon count at sea level is about 1 hit per cm2 per minute for a detector that is placed horizontally. 

The CosmicWatch Physic Paper is an essential read is you are considering building a Muon detector. Apart from describing in great detail how Muons are formed, it also outlines a number of practical experiments that can be performed.

Building the CosmicWatch detectors

Follow the build instructions (downloadable from the UKRAA and the CosmicWatch GitHub repository) and also the excellent YouTube videos listed in the instruction manual.  

Construction is very straightforward although all components are SMD units. The UKRAA sells basic kits, partially completed kits or fully assembled units. Particular care needs to be taken when covering the Scintillator block with initially reflective foil and then 5 layers (in my case) of black insulating tape. I cannot emphasise enough how important that the Scintillator and SiPM sensor face are totally light tight. The Scintillator block will 'flash' with a single photon if a charged particle passes through it. This photon is then reflected off the reflective aluminium foil lining until it lands on the SiPM sensor and triggers a count. Any stray light will at best trigger counts or at worst, overload the SiPM.

Following construction and a few weeks of 'bedding in' to make sure there were no problems, I customised the software and also 3D printed a housing to keep everything stable.

The computer is a Raspberry Pi model B (very early model and more than adequate for collecting data) and all scripts - both on the Pi and the web/data server are written in Python3.

I also made  some minor changes to the Arduino code that sits on the detectors to 'lock' one of the devices as a Slave. I have no intention of collecting data from the Master (at this time).

Changes to the Arduino Code:

The code is downloadable from the GitHub repository linked above. You need to select whether you want to use a local (on the detector) memory card or a remote computer. This determines whether you install SDCard.ino or OLED.ino as the control program. However, firstly you need to edit naming.ino first with a suitable name for your host and then upload. I wanted to store all my data on the attached Raspberry Pi so installed a standard version of OLED.ino on the Master, named MUON2, and the modified Slave-OLED.ino on MUON1. MUON1 was was installed inverted (to increase the effective capture area) and the code edited so that the display reads the rights way up.

Normally, when the units are started from scratch, the unit that powers up later (within 2 seconds) is determined to be the Slave unit. This can be very hit or miss when restarting the devices and I only wanted a single USB data cable linking to the Slave unit and I did not want to record data from the Master.

Changes to OLED.ino (Program Subset)
void setup() {
  analogReference (EXTERNAL);
  ADCSRA &= ~(bit (ADPS0) | bit (ADPS1) | bit (ADPS2));  // clear prescaler bits
  ADCSRA |= bit (ADPS0) | bit (ADPS1);                   // Set prescaler to 8
  Serial.begin(9600);
  
  display.begin(SSD1306_SWITCHCAPVCC, 0x3C);                               
  pinMode(3, OUTPUT);
  pinMode(6, INPUT);
  if (digitalRead(6) == HIGH) {
      SLAVE = 1;
      MASTER = 0;
      digitalWrite(3,HIGH);
      delay(1000);}

  else{
      delay(10);
      //MK edit to force slave mode
      MASTER = 0;
      SLAVE = 1;
      pinMode(3, OUTPUT); // was pin 6
      digitalWrite(3, HIGH);} // was pin 6, high

  if (OLED == 1){
      display.setRotation(0);         // Upside down screen (0 is right-side-up) MK - Was 2
      OpeningScreen();                // Run the splash screen on start-up
      delay(2000);                    // Delay some time to show the logo, and keep the Pin6 HIGH for coincidence
      display.setTextSize(1);}

  else {delay(2000);}
  digitalWrite(3,LOW);
  if (MASTER == 1) {digitalWrite(6, LOW);}

  Serial.println(F("#")); // to save space, the multiple # symbols have been removed. Not used in my data collection routine
  Serial.println(F("CosmicWatch:, The Desktop Muon Detector"));
  //Serial.println(F("Questions? This email address is being protected from spambots. You need JavaScript enabled to view it."));
  Serial.println(F("#"));
  Serial.println(F("Y,M,D,H,M,S,Event, Ardn_time[ms], ADC[0-1023], SiPM[mV], Deadtime[ms], Temp[C], Name")); //This header line has been altered as I edited the Python Code on the Raspberry Pi
  Serial.println(F("#"));

  get_detector_name(detector_name);
  Serial.println(detector_name);
  get_time();
  delay(900);
  start_time = millis();
  
  Timer1.initialize(TIMER_INTERVAL);             // Initialise timer 1
  Timer1.attachInterrupt(timerIsr);              // attach the ISR routine
  
}

Now, if both units are reset, the Slave unit will always be a Slave. 

One other minor change that was made to the OLED.ino code was to add a ',' to comma separate the fields in the datablock that is passed to the PC.

     if (MASTER == 1) {
          analogWrite(3, LED_BRIGHTNESS);
          sipm_voltage = get_sipm_voltage(adc);
          last_sipm_voltage = sipm_voltage; 
          Serial.println((String)count + ", " + time_stamp+ ", " + adc+ ", " + sipm_voltage+ ", " + measurement_deadtime+ ", " + temperatureC);} // added commas
  
      if (SLAVE == 1) {
          if (keep_pulse == 1) {   
              analogWrite(3, LED_BRIGHTNESS);
              sipm_voltage = get_sipm_voltage(adc);
              last_sipm_voltage = sipm_voltage; 
              Serial.println((String)count + ", " + time_stamp+ ", " + adc+ ", " + sipm_voltage + ", " + measurement_deadtime+ ", " + temperatureC);}} //added commas

 

PC code - Changes to Import_Data.py

I made extensive modifications to the Import_Data program that can be downloaded from the GitHub repository. There are a few changes that I wanted:

  • I wanted to be able to automate the program so that it would auto start without any input from the user
  • I wanted to be able to copy and input the data into a spreadsheet if needed. This meant changing the date field so that the data elements are separated.
  • I also wanted to be able to integrate the data to count number of hits per minute and then display the results graphically. A great deal of the code is unchanged so only the option that has been edited is the initialisation routine and the auto-answers set in Option 1 of the program.

Note that this will not work if two units are plugged into the Raspberry Pi using USB data cables. Only the Slave is connected using USB.

These alterations comprised of changes to the Import_Data python code and also a new Python program that sits on the Web/Data Server.

I am not going to publish the changes to the MIT Import_Data code here. Please email me if you want a copy of the code and an explanation. This email address is being protected from spambots. You need JavaScript enabled to view it.

The code that sits on the server performs an SFTP connection to the Raspberry Pi, downloads and saves the data as a CSV file (daily) and generates a graph. It's not the best code in the world, but it does appear to work.

import sys
#print(sys.path)
import pysftp
import datetime
import numpy as np
import matplotlib.pyplot as plt
import time
import os

date = datetime.datetime.now()

yesterday = date - datetime.timedelta(1) # This is run after midnight and collects the previous days data

filedate = (yesterday.strftime("%Y") + yesterday.strftime("%m") + yesterday.strftime("%d"))
filepath = "/home/pi/data/muon1" #on the Pi
filepath2 = "/home/martyn/data/muon/" #on the local server

filename="muon1_" + filedate + ".csv" #csv

fullpath = filepath2 + filename
os.chdir(filepath2) # /data/muon
#print(fullpath)
cnopts = pysftp.CnOpts()
cnopts.hostkeys = None
with pysftp.Connection('192.168.0.xx', username='xxxx', password='xxxxxxxx') as sftp:
    with sftp.cd ('data/muon1'):
#        sftp.get(fullpath)
        sftp.get(filename)

muondata = np.genfromtxt(fullpath, dtype=(int,int,int,int,int,float,int,int,int,float,int,float,str), delimiter=",", skip_header=1)
#muondata = np.loadtxt(fullpath,delimiter=",",dtype=str,skiprows=1)
#initialise the array
# array size is 1440 (min in day by two colums, (time in number of minutes and count)
s = (1440,2)
arr = np.zeros((s) , dtype=int)

f=0
while f < 1440: #Minutes in a Day
        arr[f,0]=f
        f+=1

#counts used in the loop
ct1 = 0
ct2 = 0

mdlen=len(muondata) #length in rows of the table

#print (mdlen)

while ct1 <  mdlen:
        list1 = (muondata[ct1])
        h1 = (list1[3])
        m1 = (list1[4])
        hm1 = ((h1*60)+m1)
#       print (ct1,hm1)
        ct1+=1
        if ct1 < mdlen:
                list2 = (muondata[ct1])
                h2 = (list2[3])
                m2 = (list2[4])
                hm2 = ((h2*60)+m2)

                if hm1 == hm2:
                        ct2+=1
                else:
#                       print (hm1,ct2)
                        arr[hm1,1] = ct2 # copy the muon count data into the table
                        ct2 = 0

muonmin=arr[:,0]
muoncnt=arr[:,1]
fig,ax = plt.subplots(sharex=True,figsize=(12,6))
charttitle = filedate + "  " + "-- Muon Count -- Oak Bank Observatory, Willaston UK"
ax.bar(muonmin,muoncnt, color='red', label = 'Count', lw=1)
ax.set_xlim(0,1440)
ax.set_ylim(0,250)
ax.set_xticks(np.arange(0,1440,60))
ax.set_yticks(np.arange(0,250,50))
ax.set_title(charttitle) 

# draw and save the graph
plt.ion()
plt.xlabel("Minutes (60 second data integration)")
plt.ylabel("Count")
plt.grid(True, which='both')
#plt.show()
image = "Muon"+filedate+".png"
time.sleep(2)
plt.savefig(image)
plt.close('all')


 

Finally, I 3D printed an enclosure to keep the Muons and Raspberry Pi stable to hopefully ensure reliability.