Resampling very high resolution data 4Mhz to 10 Hz
rutambhara.joshi
Hi all,
I am working with an instrument that produces very high-resolution data roughly 4 Mhz however, it is not perfect, i.e. doesn't have data point every 0.25 microseconds, that part is a bit random. I have a 2d matrix with particle counts and a time series at this time resolution.
So, I wrote a code that can average data between a time gap of 0.1 seconds. However, when I try to process 1 day of data it takes my code more than 30 minutes, is there a better way of coding this type of problem?
Many thanks,
Function pops_10Hz_grabber(pops_time, size_2d, time_grid)
wave pops_time, size_2d,time_grid // time grid is a perfact 10hz timeseries that I use to average high resolution data of pops_time and size_2d.
make/o/d/n=(numpnts(pops_time)) POPSTStore=nan
make/o/d/n=(numpnts(pops_time),11) POPSTStore_size=nan
variable i,k,j
variable cnt
variable POPSStart, POPSEnd
j=0
for(i=0;i<numpnts(time_grid)-1;i+=1)
variable chckVal = 0
variable startT= pops_time[i] // here diff betweeen startT and endT is 0.1 sec always.
variable endT = pops_time[i+1]
//print i
k = j
do
if(pops_time[k] >= startT && chckVal == 0 )
POPSStart=k
chckVal = 1 // once it finds the match chacval is 1 so the if statement doesn't fullfill and moves on to the next point
//print POPSStart
elseIf(chckVal == 1 && pops_time[k] > endT) // same way finds the endtime here in the high res timewave
POPSEnd =k-1 //
j=k // saving the position of where endT was found in J, so the next loop start from that position
k=numpnts(pops_time) // forcing k to be max value of the loop here and later adding +1 to it so I can come out of the loop
endif
k+=1
while(k<numpnts(pops_time))
duplicate/d/o/R=[POPSSTart, POPSEnd] pops_time, tempT // extracting data that fall into time gap of 0.1 sec.
wavestats/q tempT
// print V_npnts
POPSTStore[i] = V_avg
duplicate/d/o/R=[POPSSTart, POPSEnd][] POPSTStore_size tempP // not sure how to extract 2d data but I need to that too.
//matrixop/o test_data=POPSTStore
endFor
End
wave pops_time, size_2d,time_grid // time grid is a perfact 10hz timeseries that I use to average high resolution data of pops_time and size_2d.
make/o/d/n=(numpnts(pops_time)) POPSTStore=nan
make/o/d/n=(numpnts(pops_time),11) POPSTStore_size=nan
variable i,k,j
variable cnt
variable POPSStart, POPSEnd
j=0
for(i=0;i<numpnts(time_grid)-1;i+=1)
variable chckVal = 0
variable startT= pops_time[i] // here diff betweeen startT and endT is 0.1 sec always.
variable endT = pops_time[i+1]
//print i
k = j
do
if(pops_time[k] >= startT && chckVal == 0 )
POPSStart=k
chckVal = 1 // once it finds the match chacval is 1 so the if statement doesn't fullfill and moves on to the next point
//print POPSStart
elseIf(chckVal == 1 && pops_time[k] > endT) // same way finds the endtime here in the high res timewave
POPSEnd =k-1 //
j=k // saving the position of where endT was found in J, so the next loop start from that position
k=numpnts(pops_time) // forcing k to be max value of the loop here and later adding +1 to it so I can come out of the loop
endif
k+=1
while(k<numpnts(pops_time))
duplicate/d/o/R=[POPSSTart, POPSEnd] pops_time, tempT // extracting data that fall into time gap of 0.1 sec.
wavestats/q tempT
// print V_npnts
POPSTStore[i] = V_avg
duplicate/d/o/R=[POPSSTart, POPSEnd][] POPSTStore_size tempP // not sure how to extract 2d data but I need to that too.
//matrixop/o test_data=POPSTStore
endFor
End
Joshi
Hello Joshi,
I am a bit confused about your data storage. Could you please explain what you are storing in the 2D wave?
AG
June 14, 2023 at 11:57 am - Permalink
Hi AG,
So, that's the part I haven't figured out yet. I have a 2d wave with about 4Mhz resolution, and I am trying to average that 2D so it reduces to 10 hz. I have 11 columns in my 2d wave and like to apply the same type of averaging to each column and store this 2D wave in a reduced-size matrix.
Many thanks,
Joshi
June 14, 2023 at 01:26 pm - Permalink
You can run your loop over the time wave to determine the POPSSTart, and POPSEnd rows that are relevant for each time interval. You can then use MatrixOP subRange(w,POPSSTart,POPSEnd,0,maxCol) to pick the relevant rows. Finally, you wrap this up with MatrixOP averageCols(). This would be a one-line MatrixOP command in your loop that will look like:
MatrixOP/O POPSTStore_size=setRow(POPSTStore_size,i,averageCols(subRange(w,POPSSTart,POPSEnd,0,maxCol)))
I hope this helps,
A.G.
June 14, 2023 at 05:11 pm - Permalink