Welcome, Guest

Want to take part in these discussions? Sign in if you have an account, or apply for one below

• Forgot your password?
• Apply for membership

1. • CommentAuthorsyrex314
   • CommentTimeApr 30th 2010

    
   I'm trying to calculate a partial cross-correlation, looping through two arrays simultaneously. I find that if I access the elements by ILArray.GetValue(i), it is far slower than first copying the whole array to a double[], and then accessing the individual elements.
   
   E.g. the following code is slow:
   ILArray<double> x;
   ILArray<double> y;
   ...
   for( int i = 0; i < x->Length; i++ )
   for( int j = 0; j < x->Length; j++ )
   ...
   something = x.GetValue(i) * y.GetValue(i+j);
   
   
   and the following code is much faster:
   
   ILArray<double> x;
   ILArray<double> y;
   double[] xx = new double[x->Length];
   double[] yy = new double[y->Length];
   x.ExportValues( ref xx );
   y.ExportValues( ref yy );
   ...
   for( int i = 0; i < x->Length; i++ )
   for( int j = 0; j < x->Length; j++ )
   ...
   something = xx[i] * yy[i+j];
   
   
   Am I going about this the wrong way? It seems to me that even in the second example, things could be faster if I didn't have to copy each ILArray into a double[]. I know this violates principles of abstraction, but is there any way at all to gain access to the m_data internal storage of ILArray?
   
   Thanks.
2. • CommentAuthorsyrex314
   • CommentTimeApr 30th 2010

    
   I think I found my answer: ILArray.InternalArray4Experts[].
3. • CommentAuthorhaymo
   • CommentTimeApr 30th 2010

    

   hm, like the name says: this array is to be used with extreme care! It is by no means garanteed, that it reflects the underlying elements in the way represented by the ILArray itself !!!! For example, if you build a reference array from ILArray<T> A by - lets say- transposing A, than A.T is a reference array. the ordering of A.T elements is different for matrices than that of A - even if both share the same underlying  "InternalArray4Experts" System.T[] array.

   In general, a more better way is to work around adressing single elements at all. There will most probably be a way to get rid of any loops and work on the ILArrays itself. ("A * B"). Post your code, if you need any help on that.

   PS: One more note on "InternalArray4Experts": technically spoken, there is nothing wrong by accessing this System.Array - as long as the following restrictions are hold at least:

   * make sure to only read from that underlying data array! There might be other ILArrays using this storage. The result of writing to it would be corruption of those arrays.

   * make sure, the "InternalArray4Experts" reflects the storage of the ILArray correctly. This can be done by "ILArray.Detach()" which will change any referencing array into an solid array.

   * anyway - trying not to use the "InternalArray4Experts" will always be safer :)

4. • CommentAuthorsyrex314
   • CommentTimeApr 30th 2010

    
   Below I have pasted the full code for my function. I have data coming in to inArray, and I want to use an autocorrelation to calculate how many elements it needs to be shifted by in order to align the data to an accumulation buffer for averaging. If the incoming data is not aligned to the accumulation buffer, then high-frequency features will get washed out in the average.
   
   The reason I want to iterate over indices instead of simply calculating the vector inner product of inArray and accumulationBuffer, is that I I know the two are going to be closely aligned. When I implemented this in C++, I found the only way to get this to keep up with the data stream I am working with is to restrict the cross-correlation shift range from -n/16 to +n/16. If there's a better way to approach this, I'm all ears!
   
   In the meantime, I'll take your advice and call Detach() at the top, and make sure to not write anything to the underlying System.Array.
   
   private int FindCorrelationShift( ref ILArray<double> inArray ) {
   // calculates how many elements to shift a data array by to line up with the accumulation array
   // positive index shifts to the left, discarding earlier indices
   
   // first find the mean of both arrays
   double meanReference = 0;
   double meanData = 0;
   
   int i;
   int n = accumulationBuffer.Length;
   
   meanReference = (double)ILMath.mean(accumulationBuffer);
   meanData = (double)ILMath.mean(inArray);
   
   // compute the partial cross-correlation
   double ccmax = 0; // peak value of cross correlation
   double cc; // value of cross-correlation at current shift
   int iShift = 0; // current shift index
   for( int j = -n/16; j < n/16; j++ ) {
   // shift data by i
   cc = 0;
   
   for( i = 0; i < n/4; i++ ) {
   // not the most efficient way to do this, but it should work to just check bounds first
   // the better way is to think about indices carefully
   if (((i + j) >= 0) && ((i + j) < n))
   cc += (accumulationBuffer.InternalArray4Experts[i] - meanReference)
   * (inArray.InternalArray4Experts[j + i] - meanData);
   }
   
   // check for correlation peak value
   if( cc > ccmax ){
   ccmax = cc;
   iShift = j;
   }
   }
   
   return iShift;
   
   }
5. • CommentAuthorhaymo
   • CommentTimeMay 1st 2010

    

   I see... Dont know exactly, why and how this restriction you are using comes into play. But the "true" cross corellation is much faster implemented via fft. I tested your code against the fft version and the speed up seems to be significant.

   using System;
   using System.Collections.Generic;
   using System.Text;
   using ILNumerics;
   using ILNumerics.BuiltInFunctions;
   
   namespace ConsoleApplication1 {
       class Program {
           static void Main(string[] args) {
   
               int len = 50000;
               ILArray<double> accBuf = Computation.createTestData(len);
               int shift = 102;
               ILArray<double> data = accBuf["0;" + shift + ":end,0:" + (shift - 1)];
               int sFFT = Computation.FindCorrelationShiftFFT(data, accBuf);
               int s = Computation.FindCorrelationShift(data, accBuf);
               System.Diagnostics.Debug.Print(String.Format(
                   "OrigShift:{0} FoundShift:{1} FoundViaFFT:{2}", shift, s, sFFT));
           }
   
           private class Computation : ILMath {
   
               public static int FindCorrelationShift(
                   ILArray<double> inArray, ILArray<double> accumulationBuffer) {
                   // calculates how many elements to shift a data array by to line up with the accumulation array
                   // positive index shifts to the left, discarding earlier indices
                   // ... THIS IS YOUR CODE ABOVE ...

                   return iShift;
   
               }
   
               public static int FindCorrelationShiftFFT(
                   ILArray<double> A, ILArray<double> acc) {
                  
                   ILArray<double> cc = ifftsym(conj(fft(A)) * fft(acc));
                   ILArray<double> shift = empty();
                   max(cc, ref shift, A.Dimensions.FirstNonSingleton());
                   return (int)shift.GetValue(0);
               }
   
               public static ILArray<double> createTestData(int len) {
                   return rand(1, len) + sin(linspace(0, pi, len));
               }
           }
       }
   }

   See the full console application file here.

   If you need to "restrict" the search for the shift to some fraction of the whole length, this could easily be done in the max(..) statement: cc["0:"+limit+"," + endlimit + ":end"] (or similar, not tested!)

6. • CommentAuthorsyrex314
   • CommentTimeMay 1st 2010

    
   Thanks! I'll give the FFT version a try. I'm embarrassed to admit I forgot about this approach from my undergrad DSP class.
   
   By the way, I'm very impressed with ILNumerics. I've been looking for a way to break free from Matlab for a while now.
7. • CommentAuthorhaymo
   • CommentTimeMay 1st 2010

    

   Thank you. One more note to the "InternalArray4Experts". Whenever such problems come up, one should try to get around it in the following order: 

   1) "Unloop" the code - prevent loops over any elements. This can most often be done be reformulating the algorithm.

   2) Check, if you can utilize one of the iterators. The ILArray<T>.Values collection iterates over all elements in a safe manner and gives read-only access to them. Also, it is faster than querying elements via GetValue()

   3) When single element access is needed nevertheless:

   3.1) Detach the array at first. This may also speed up the GetValue() functions.

   3.2) Use GetValue(), specify only one parameter if possible, (ie. sequential access on one dimension)

   3.3) Export the array into a predefined System.Array and work on that

   3.4) Use "InternalArray4Experts" with lots of care...

    

8. • CommentAuthorsyrex314
   • CommentTimeJun 7th 2010

    
   Haymo,
   
   Thanks for your help. I finally got around to implementing the FFT cross-correlation in my application, and it works very well! I also had to implement a circshift() function to line up the arrays once the cross-correlation calculates the shift. I've pasted it here if anyone may find it to be useful:
   
   public static ILArray<double> circshift(ILArray<double> inArray, int nShift)
   {
   //circularly shift an array to the right.
   ILArray<double> outArray;
   
   if (nShift < 0)
   {
   // turn a left shift into a right shift of n - nShift
   nShift = inArray.Length + nShift;
   }
   
   
   if (nShift == 0)
   outArray = inArray;
   else
   {
   outArray = ILMath.zeros(inArray.Length, 1);
   
   // make the first nShift elements of outArray equal to the last nShift elements of inArray
   outArray["0:" + Convert.ToString(nShift - 1)]
   = inArray[Convert.ToString(inArray.Length - nShift)
   + ":end"];
   
   // make the last (n-nShift) elements of outArray equal to the first (n-nShift) of outArray
   outArray[Convert.ToString(nShift) + ":end"]
   = inArray["0:" + Convert.ToString(inArray.Length - nShift - 1)];
   }
   return outArray;
   
   }
9. • CommentAuthorhaymo
   • CommentTimeJun 7th 2010

    

   thanks :)

   just one more note: be carefull on assignments like those:

   outArray = inArray,

   since both will point to the same object afterwards! If this is not required behaviour, use

   outArray = inArray.C  to clone the array or outArray = inArray.R  to create a reference instead.