SignalMixin.smooth_frames: assorted optimizations.

The variable tspf is not used and can be removed.

sig_len can be calculated using // rather than / and truncating to an
integer (the latter could be incorrect for records longer than 2**53.)

The output array can be allocated using np.empty instead of np.zeros
since it will be fully initialized by the subsequent loop.

Since each frame is independent of the others, the loop can be broken
up into blocks (here, arbitrarily chosen as 2**16 frames) to reduce
temporary memory usage while still spending most of our CPU time in
numpy operations.

Author: Benjamin Moody

wfdb/io/_signal.py

@@ -827,9 +827,6 @@ def smooth_frames(self, sigtype='physical'):
             if spf[ch] is None:
                 spf[ch] = 1
 
-        # Total samples per frame
-        tspf = sum(spf)
-
         # The output data type should be the smallest type that can
         # represent any input sample value.  The intermediate data type
         # must be able to represent the sum of spf[ch] sample values.
@@ -854,7 +851,7 @@ def smooth_frames(self, sigtype='physical'):
             raise ValueError("sigtype must be 'physical' or 'digital'")
 
         n_sig = len(expanded_signal)
-        sig_len = int(len(expanded_signal[0])/spf[0])
+        sig_len = len(expanded_signal[0]) // spf[0]
         input_dtypes = set()
         for ch in range(n_sig):
             if len(expanded_signal[ch]) != sig_len * spf[ch]:
@@ -868,15 +865,22 @@ def smooth_frames(self, sigtype='physical'):
             if all(dt <= output_dtype for dt in input_dtypes):
                 break
 
-        signal = np.zeros((sig_len, n_sig), dtype=output_dtype)
+        signal = np.empty((sig_len, n_sig), dtype=output_dtype)
+
+        # Large input arrays will be processed in chunks to avoid the need
+        # to allocate a single huge temporary array.
+        CHUNK_SIZE = 65536
 
         for ch in range(n_sig):
             if spf[ch] == 1:
                 signal[:, ch] = expanded_signal[ch]
             else:
                 frames = expanded_signal[ch].reshape(-1, spf[ch])
-                signal_sum = np.sum(frames, axis=1, dtype=intermediate_dtype)
-                signal[:, ch] = signal_sum / spf[ch]
+                for chunk_start in range(0, sig_len, CHUNK_SIZE):
+                    chunk_end = chunk_start + CHUNK_SIZE
+                    signal_sum = np.sum(frames[chunk_start:chunk_end, :],
+                                        axis=1, dtype=intermediate_dtype)
+                    signal[chunk_start:chunk_end, ch] = signal_sum / spf[ch]
 
         return signal