I am quite new to Mathematica and I am looking for some advice in order to improve and optimize my code.

Here is what I am trying to do:

  1. Open a stream

  2. Read the stream line by line and store it in "data"

  3. For each value of "data", I would like to find the position of the first peak and store it in "firstPeaks"

  4. Close the stream

How could I optimize this process and make the findAllPeaks function run faster? I tried including the parallelization option but that didn't really change the run time. Running the code as is takes about 15 seconds for a $ 2500 \times 3000 $ array (2500 being the number of data points and 3000 being the number of bytes per data point), which is significantly smaller than what I would like to work with.

I did read the following article: Jon McLoone 10 Tips for Writing Fast Mathematica Code and I did try to follow some of the suggestions, namely:

  • Using compile.
  • Using built-in functions such as Apply, PeakDetect and Position.
  • Using parallelization, which didn't seem to help.
  • Using Block, which also didn't seem to help much in this case.
  • And using "EvaluateSymbolically" -> False to avoid any symbolic calculations.
readlines[n_, m_] := Block[{}, SetStreamPosition[str, index[[n]]]; ReadList[str, {Byte}, m]]

findPeaks = Compile[{{d, _Integer, 1}}, PeakDetect[d, .1, .01], RuntimeOptions -> {"EvaluateSymbolically" -> False}]

str = OpenRead[FileNameJoin[{path, fileName}], BinaryFormat -> True];

index = Table[pos = StreamPosition[str]; Skip[str, Byte, lineLength];
   pos, {streamLength}];

Do[data[[i]] = Flatten[readlines[i, lineLength]], {i, 1, streamLength}];

Do[firstPeaks[[i]] = FirstPosition[Apply[findAllPeaks, data[[i]]], 1., {0}], {i, 1, streamLength}]; // AbsoluteTiming

  • 1
    $\begingroup$ Did you try using Import instead of low-level file reading functions? Also, PeakDetect is not in the list of compilable functions which means that putting it inside Compile will not do anything (it may even slow things down a little.) $\endgroup$ – C. E. Apr 7 '19 at 6:12
  • 1
    $\begingroup$ @etotheix Welcome to Mathematica.StackExchange. I notice that your question is well written. However, it would be easier to answer if you could provide also an example file or an example data set. Apart from that, I second C.E.'s suggestions: try to read the data in one chunk as a matrix, e.g. with Import. Afterwards Map or ParallelMap PeakDetect over it. $\endgroup$ – Henrik Schumacher Apr 7 '19 at 8:17
  • 1
    $\begingroup$ I suspect that PeakDetect uses GaussianFilter under the hood; the documentation indicates that. For some reason that is not clear to me, GaussianFilter performs poorly without setting the option Method -> "Gaussian". Unfortunately, PeakDetect does not feature such an option. In this post, I showed how GaussianFilter can be sped up with ListConvolve. One could build a PeakDetect from scratch with the information there. It is only a matter of guessing the parameters and the settings for Padding etc. correctly. $\endgroup$ – Henrik Schumacher Apr 7 '19 at 8:56
  • 4
    $\begingroup$ I implemented PeakDetect from scratch in my answer here. $\endgroup$ – C. E. Apr 7 '19 at 9:16
  • 2
    $\begingroup$ @HenrikSchumacher Does SetOptions[GaussianFilter, Method -> "Gaussian"] have the desired effect on PeakDetect? (Trace[ PeakDetect[data], _GaussianFilter ] does confirm GaussianFilter is called. $\endgroup$ – Michael E2 Apr 7 '19 at 13:55

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.