US20160089048A1 - Time transformation of local activation times - Google Patents

Abstract

An automatic method of determining local activation time (LAT) from at least four multi-channel cardiac electrogram signals including a ventricular channel, a mapping channel and a plurality of reference channels. The method comprises (a) storing the cardiac channel signals, (b) using the ventricular and mapping channel signals and a first reference channel signal to compute LAT values at a plurality of mapping-channel locations, (c) monitoring the timing stability of the first reference channel signal, and (d) if the timing stability of the monitored signal falls below a stability standard, using the signal of a second reference channel to determine LAT values. Substantial loss of LAT values is avoided in spite of loss of timing stability.

Description

FIELD OF THE INVENTION
• This invention is related generally to the field of electrophysiology, and more particularly to technology for accurate measurement of parameters within body-surface ECG, intracardiac, and epicardial electrical signals such as heart rates and local activation times and the assessment of the quality of such measurements.
BACKGROUND OF THE INVENTION
• The invention disclosed herein involves the processing of multiple channels of electrical signals which are produced by the heart. These channel signals include the ECG signals from body-surface electrodes and signals from electrodes within the body, i.e., intracardiac signals from within vessels and chambers of the heart and epicardial signals from the outer surface of the heart. Throughout this document, the term “multi-channel cardiac electrogram” (or “MCCE”) is used to refer to all of these types of channels, and when specific types are appropriate, specific nomenclature is used. This new terminology (MCCE) is used herein since the term “ECG” sometimes only refers to body-surface measurements of cardiac performance.
• A major component in cardiac interventional procedures such as cardiac ablation is the display of cardiac data which is extracted from the MCCE signals captured by an array of electrodes placed on the body surface and within and on the structures of the heart itself. Among the important data which are displayed are ventricular pulse interval (time between heart beats), intracardiac cycle length (time between the activations in arrhythmias (such as atrial fibrillation), relative time differences between related activations in two intracardiac channels to generate activation maps, and assessments of signal strength, variability, and other measures of signal quality within MCCE signals.
• Cardiac interventional electrophysiology procedures (e.g., ablation) can be extremely time-consuming, and the reliable determination and presentation of such cardiac parameters is an important element in both the quality of the procedures and the speed with which they can be carried out. Often the data which are presented to the electrophysiology doctor during such procedures exhibit high variability contributed not only by the performance of the heart itself but by unreliable detection of certain features of the MCCE signals. Therefore there is a need for more reliable and more rapid algorithms to process body surface and intracardiac signals obtained during an electrophysiology (EP) procedure.
• MCCE electrodes capture the electrical signals in the cardiac muscle cells. As mentioned above, some MCCE electrodes may be attached to the body surface (ECG) and some may be positioned inside cardiac veins, arteries and chambers (intracardiac) and on the outer surface of the heart (epicardial) as conductive elements at the tips or along the lengths of catheters introduced into the body and maneuvered into position by the EP doctor. The electrical signals within the heart muscles and which flow therefrom to other regions of the body have very low voltage amplitudes and therefore are susceptible to both external signal noise and internally-generated electrical variations (non-cardiac activity). In addition, cardiac arrhythmias themselves may be highly variable, which can make reliable extraction of cardiac parameters from MCCE signals difficult.
• One important cardiac parameter used during such procedures is the time difference between the activations occurring within two channels, both of which contain the electrical signals of an arrhythmia. This measurement is called local activation time (LAT), and measurement of a plurality of values of LAT is the basis for the generation of an activation map. The map displays information about the sequence of activations of cardiac muscle cells relative to each other, and this sequence of information is combined with physical anatomical position information to form the map. An activation map then provides guidance to the EP doctor for the process of applying therapies to heart muscle cells which can terminate cardiac arrhythmias and permanently affect the heart to prevent recurrence of such arrhythmias.
• The entire process of determining LAT is referred to as mapping because all of the information generated by analysis of the MCCE signals is combined in a single computer display of a three-dimensional figure that has the shape of the heart chamber of interest and display of additional image qualities such as color that convey the sequence of electrical activity (activation map) or possibly other qualities of the electrical activity (e.g., voltage map). These images are similar in style to weather maps common today in weather-forecasting. Such a cardiac map becomes a focus of attention for the EP doctor as he directs the motion of catheters in the heart to new positions, and an algorithm which processes the MCCE signals produces measurements from the electrodes in new positions. As this process continues, the map is updated with new colored points to represent additional information about the electrical activity of the heart.
• During a mapping procedure, the timing relationships of muscle depolarizations typically must be determined for hundreds of locations around a heart chamber which may be experiencing an abnormal rhythm. The locations are often examined, one at a time, by moving an exploring cardiac-catheter electrode (mapping-channel electrode) from location to location, acquiring perhaps only a few seconds of signal data at each location. To compare timing relationships, a different electrode (reference-channel electrode) remains stationary at a single location and continuously acquires a reference signal of the rhythm.
• The collection of timing relationships and anatomical locations constitutes an activation map (LAT map). As described above, a relatively large number of individual LAT values are used to generate a useful LAT map. Many different locations may serve adequately as alternative reference locations, but it has been critical in the present state-of-the-art that whatever location is used as the reference, one activation map is committed for the entire duration to that reference location only.
• U.S. patent application Ser. No. 13/922,953 (Brodnick) filed on Jun. 20, 2013 discloses several aspects of improved methods for determining LAT. (Such application and the invention of the present application are commonly-owned, and Donald Brodnick is also an inventor of the present invention.) The Brodnick application discloses LAT-determination methods which include replacement of cardiac channels when the quality of such channel signals falls below a standard measure of channel-signal quality. Much of the disclosure of the Brodnick application is included herein since it provides excellent background information for the improved LAT-determination methods disclosed herein.
• Occasionally a reference electrode is bumped or becomes disconnected. In these cases, additional data cannot be collected to extend the map (add more LAT values to the map) because the timing relationships are no longer comparable (based on the same reference channel signal). The doctor either makes his interpretation of the map based on an incomplete map or establishes a new reference and begins to create a new map, having lost the time and effort which to this point in the procedure had been expended. At a few seconds of signal acquisition per location, a few seconds of catheter motion between locations, and hundreds of locations, the amount of time and effort wasted if a map must be restarted can be very significant. Furthermore, extending the total procedure time adds more risk of complications for the patient.
• Because the heart is constantly contracting and other catheters are continually being repositioned, a procedure may last for several hours, during which time the patient even may need to be moved. Occasionally a reference electrode either makes poor contact or may shift position, in which case the constant timing relationship is disrupted (timing stability is lost) and additional locations cannot be studied in relationship to accumulated data. As described above, the resulting incomplete activation map may be worthless, requiring a new map, extending the procedure and adding cost and risk to the patient.
• Thus there is a need for an automatic method of determining local activation time (LAT) from multi-channel cardiac electrogram signals which avoids substantial loss of LAT values in spite of losses of timing stability in reference channels during a local activation time mapping procedure.
• The generation of position information and its combination with cardiac timing information is outside the scope of the present invention. The focus of the present invention is the processing of MCCE signals to measure time relationships within the signals, the two most important of which are cycle length (CL) and local activation time (LAT).
• Currently-available MCCE-processing algorithms are simplistic and often provide inaccurate measurements which cause the activation map and many other cardiac parameter values to be misleading. A misleading map may either (1) compel the doctor to continue mapping new points until apparent inconsistencies of the map are corrected by a preponderance of new, more-accurately measured map points or (2) convince the doctor to apply a therapy to a muscle region which actually makes little or no progress in the termination of an arrhythmia, again prolonging the procedure while the EP doctor maps more points in an attempt to locate new regions where therapy may be effective.
• Currently, computer systems which assist doctors in the mapping process have manual overrides to allow a technician, or sometimes the EP doctor himself, to correct the measurements made by the system automatically. This requires a person to observe a computer display called the “Annotation Window” which shows a short length of the patient's heart rhythm, perhaps 3-5 heart beats as recorded in 3-8 channels (signals from MCCE electrodes).
• The channels of the annotation window are of several types. Usually there will be a body surface ECG lead such as lead II that identifies when ventricular activity is occurring. (It is also possible that the ventricular activity may be sensed by an intracardiac channel electrode.) There is one channel, identified as a reference channel, the electrode of which ideally remains in a fixed position during the entire map-generating procedure. There is at least one other intracardiac channel (the mapping channel) which senses the electrical signal at a catheter tip, the precise three-dimensional position of which is determined by other means. The electrical activity in the mapping channel is compared to the activity in the reference channel to determine the local activation time (LAT) which is used to color the map at that precise three-dimensional position.
• Intracardiac channels may be of either the bipolar or unipolar recording types, and the inventive measurement method disclosed herein can be applied to both types of signals. Also, since it is possible during arrhythmias for some chambers of the heart to be beating in a rhythm different from other chambers of the heart, the annotation window often contains additional channels to aid the doctor's interpretation of the data presented.
OBJECTS OF THE INVENTION
• It is an object of this invention, in the field of electrophysiology, to provide an automatic method for accurate measurement of several parameters which characterize MCCE signals.
• Another object of this invention is to provide an automatic method for such measurements which operates rapidly enough to not hinder an electrophysiologist performing procedures which utilize such a method.
• Another object of this invention is to provide an automatic method for rapid and reliable measurement of cardiac parameters to reduce the length of time certain cardiac procedures require and also reduce the X-ray exposure times for the patients.
• Another object of this invention is to provide an automatic method for rapid and reliable measurement of local activation times which are provided for the rapid generation of local activation time maps, determining the precise phase relationship between a reference channel and a mapping channel.
• Still another object of this invention is to provide an automatic method for cardiac parameter measurement which can be used in real time during certain interventional cardiac procedures.
• Another object of this invention is to provide an automatic method for rapid and reliable activation mapping which can continue providing LAT measurement when a reference signal degrades in timing stability such that it is no longer useable as a reference signal.
• Another object of the invention is to provide an automatic method for measuring cardiac parameters which is largely insensitive to the amplitude of the MCCE signals and almost entirely dependent on the timing information contained in such signals.
• Another object of this invention is to provide an automatic method for measurement of local activation times which avoids the loss of LAT values which have been determined prior to a loss of timing stability of the reference channel signal used to determine such LAT values.
• Yet another object of the present invention is to provide an automatic method for generating a single map throughout an LAT mapping procedure even when all of a plurality of reference channel signals fail intermittently at different times as long as at least one reference channel is functioning properly at any time during the mapping procedure.
• And yet another object of the inventive method is to provide reliable and accurate automatic determination of cardiac-channel timing stability and signal quality.
• These and other objects of the invention will be apparent from the following descriptions and from the drawings.
SUMMARY OF THE INVENTION
• The term “digitized signal” as used herein refers to a stream of digital numeric values at discrete points in time. For example, an analog voltage signal of an MCCE channel is digitized every millisecond (msec) using an analog-to-digital (A/D) converter to generate a series of sequential digital numeric values one millisecond apart. The examples presented herein use this sampling rate of 1 kHz, producing streams of digital values one millisecond apart. This sampling rate is not intended to be limiting; other sampling rates may be used.
• The term “velocity” as used herein refers to a signal the values of which are generally proportional to the time-rate-of-change of another signal.
• The term “velocity-dependent signal” as used herein refers to a set of possible signals which relate to the velocity of a channel signal, and in particular, retain certain properties of channel velocity. Channel signals are filtered to generate velocity-dependent signals which contain signal information which does not lose either the positive or negative activity in a channel signal. One such velocity-dependent signal is the absolute value of channel velocity; such a velocity-dependent signal is used in some embodiments of the inventive method to preserve the magnitude of the activity in a signal. Other possible velocity-dependent signals are even powers of velocity (squared, 4^thpower, etc.) which retain both the positive and negative signal activity in a velocity signal—the relative magnitudes are not critical in the present invention as long as both positive and negative activity within the signals are not masked by the filtering. Numerous other possible filtering strategies may be used to generate velocity-dependent signals, such as comparison of positive portions of the velocity with respect to a positive threshold and similarly, comparison of negative portions of the velocity with respect to a negative threshold. With respect to their use in the present invention, all velocity-dependent signals as defined herein are fully equivalent to absolute-value velocity filtering in every relevant respect.
• The term “two differenced sequential boxcar filters” as used herein refers to two boxcar filters which operate in tandem and then the difference between the two boxcar filter values is computed. Such a filtering operation is one embodiment by which a low-pass filter followed by a first-difference filter is applied. Two differenced sequential boxcar filters are illustrated inFIG. 3A and described in detail later in this document.
• The term “dot-product autocorrelation” as used herein refers to a mathematical operation applied to a time series of digital values, and this operation is generally a signal-processing application of conventional autocorrelation. Applying conventional autocorrelation to a fixed-length time series of numeric values x_igenerates another series of numeric values a_jwhich represents how well the signal x_icorrelates with itself as a function of the time difference between the signal and the same signal displaced in time by a period of time called lag. In conventional autocorrelation of a fixed-length signal x_ihaving n values in a time interval,
• 
  a_j=Σ(x_i·x_i-j)
• where the symbol Σ indicates the sum over all n−j values of x_i, represents values of time, and j represents values of lag. As used herein, the dot-product autocorrelation may be adjusted by a scale factor K as a computational convenience, in which case,
• 
  a_j=K·Σ(x_i·x_i-j)
• again where the symbol Σ indicates the sum over all n−j values of x_i, represents values of time, and j represents values of lag. Such an adjustment is not intended to be limiting to the meaning of the term. The maximum value of a_jis, of course, a₀since at lag=0, the signal perfectly correlates with itself. One such form of scale factor may include a₀such that K=k/a₀where k is a constant and its value is set for computational convenience.
• The term “magnitude-coincidence autocorrelation” as used herein refers to a modification of dot-product autocorrelation. As used herein, magnitude-coincidence autocorrelation operates on signals which first have been rectified (an absolute-value filter has been applied). Each numeric value of a fixed-length time series x_i(all values of x_i≧0) is replaced by a 1 if the value x_iis equal to or greater than a threshold value T_ACand by a 0 if x_iis less than T_AC. Further, threshold T_ACis set at some multiple of the median of all n values of x_iin the fixed-length time series. Rectified cardiac signals such as those processed by the present invention contain noise which is typically substantially smaller than the peaks within such signals. Furthermore, over all n values of such a signal, a large number of values will be close to the noise level since there are substantial periods of time between signal (electrical events) representing a heart beat. Therefore, if threshold T_AC=p·median(x_i) and p is, e.g., 4, threshold T_ACwill be just above the baseline noise in the signal x_i, and the thresholded signal X_iwill be equal to 1 only if a signal value is present which is generally not noise. Then, the magnitude-coincidence autocorrelation will have peaks for values of lag at which the time-distribution of the noise-free signal aligns (correlates well) with itself. Magnitude-coincidence autocorrelation is particularly useful when the “time” information in a signal is of more interest than the “shape” or amplitude information in a signal.
• The term “normal median” as used herein refers to the numeric value determined from a set of numeric values, such numeric value (median) being computed according to the commonly-understood mathematical meaning of the term median. The normal median of a finite set of numeric values can be determined by arranging all the numeric values from lowest value to highest value and picking the middle value from the ordered set. If there is an even number of numeric values in the set, the normal median is defined to be the mean of the two middle values of the ordered set.
• The term “set-member median” as used herein refers to the numeric value determined from a set of numeric values in a manner modified from the above-described method of median determination. In this modified determination, if there is an even number of numeric values in the set, the set-member median is either one of the two middle values in the ordered set such that the set-member median is always a member of the set of numeric values. As a practical matter, in almost all sets of real data, there are a very large numbers of data values near the median, and there is little if any difference between the two middle values.
• The term “intracardiac channel” as used herein refers to a channel of a set of MCCE signals which is connected to an internal lead, i.e., connected to a internal-surface electrode such as is at the end or along the tip of a cardiac catheter. For example, such an electrode may be in a blood vessel or in a chamber of a heart.
• The term “ventricular channel” as used herein refers to a channel of a set of MCCE signals which exhibits the dominant response of the ventricles. This may most often be a channel which is connected to an external lead, i.e., connected to a body-surface electrode. An epicardial or intracardiac channel may also sometimes be a ventricular channel.
• The term “activation” as used herein refers to a time segment within an MCCE signal which represents the passage of a depolarization wavefront within muscle cells adjacent to an MCCE electrode. An activation may sometimes be referred to as an activity trigger. Note that the terms “activations” and “activation times” may herein be used interchangeably since each activation has an activation time associated with it.
• The term “cycle length” as used herein refers to the time between neighboring activations in an MCCE signal, particularly in a reference-channel or mapping-channel signal. As used herein, the term “pulse interval” is used to connote the cycle length for a ventricular channel. The terms “ventricular pulse interval” and “intracardiac cycle length” are used to distinguish between these two measures of repetitive signals. For example, if a cardiac patient is in a period of atrial fibrillation or flutter, there may be a significant difference between the rate of occurrence of electrical events in a ventricular channel and in some intracardiac channels. The ventricular cycle length, herein called ventricular pulse interval to further distinguish it from intracardiac cycle length, may be two or three times as long as the intracardiac cycle length.
• As used herein, the terms “method” and “process” are sometimes used interchangeably, particularly in the description of the preferred embodiment as illustrated in the figures. The algorithms described as embodiments of the inventive automatic method of measuring parameters of multi-channel cardiac electrogram signals are presented as a series of method steps which together comprise processes.
• As used herein, the terms “signal” and “channel” may be used interchangeably since the inventive automatic method described herein uses signal values in the channels of MCCE signals. For example, often as used herein, the term “channel” implies the addition of the word “signal” (to produce “channel signal”) but for simplicity and textual flow, the word “channel” is used alone.
• The term “timing stability” as used herein refers to the degree to which a timing parameter, such as LAT, changes from one value to the next value during a cardiac procedure, based on a standard for timing stability. For example, an LAT may be said to be stable if has not changed from its past value (or a composite of past values) by more than a predetermined percentage or by more than a multiple of its standard deviation. Measurement of a timing parameter may of course also be affected by noise in one or more of the MCCE signals such that a determination of such parameter is degraded beyond usefulness. Such an occurrence will also be seen as a loss of timing stability.
• The term “substantial loss of LAT values is avoided” as used herein refers to largely preventing the loss of the time and effort invested by the EP doctor in capturing LAT values and not narrowly to whether or not specific numerical values for LAT are retained. Avoiding substantial loss of LAT values may mean (a) that specific LAT values are used in an unchanged form, (b) that specific LAT values are corrected in order to be useful, and/or (c) that specific LAT values are replaced by other LAT values determined from already-existing cardiac electrogram signals. In all of these situations, the LAT values, whether in changed or unchanged form, are still available to be used by the EP. Changed LAT values are herein referred to as having been transformed.
• The term “base reference channel” as used herein refers to the reference channel used in an LAT computation. LAT is computed using a mapping channel, a ventricular channel and a reference channel, and the reference channel in this group of three channels is sometimes referred to herein as the base reference channel.
• The term “signal characteristic” as used herein refers to a metric of a signal by which differences between signals may be distinguished.
• The term “center-of-power frequency” as used herein refers to the first moment of power computed from a signal frequency spectrum.
• The term “frequency-selective aggregate magnitude” as used herein refers to a value formed by combining multiple Haar transformation coefficients having differences based on the same time scale into a single value.
• The present invention is an automatic method of determining local activation time (LAT) from at least four multi-channel cardiac electrogram signals which include a ventricular channel, a mapping channel and a plurality of reference channels. The method comprises: (a) storing the cardiac channel signals; (b) using the ventricular and mapping channel signals and a first reference channel signal to compute LAT values at a plurality of mapping-channel locations; (c) monitoring the timing stability of the first reference channel signal; and if the timing stability of the monitored signal falls below a stability standard, using a second reference channel signal to determine LAT values and avoid substantial loss of LAT values in spite of loss of timing stability.
• Some preferred embodiments of the inventive automatic LAT-determining method include computing one or more timing offsets using pairs of the plurality of reference channel signals, a timing offset being LAT_K(J), the local activation time of a reference channel J based on a reference channel K and used to transform an LAT value based on reference channel J to an LAT value based on reference channel K.
• In certain preferred embodiments, using the second reference channel signal to determine LAT values includes transforming future LAT values such that they are based on the first reference channel. In some of these embodiments, LAT₂(M) is a future LAT value of mapping channel M based on the second reference channel, and future transformed values LAT₁(M) of mapping channel M based on the first reference channel are equal to a timing offset LAT₁(2) plus LAT₂(M).
• Some other preferred embodiments using the signal of a second reference channel to determine LAT values by transforming LAT values include transforming past LAT values such that they are based on the second reference channel. In some of these embodiments, LAT₁(M) is a past LAT value of mapping channel M based on the first reference channel, and past transformed values LAT₂(M) of mapping channel M based on the second reference channel are equal to a timing offset LAT₂(1) plus LAT₁(M).
• In some highly-preferred embodiments, the one or more timing offsets are computed at a plurality of times, and the value of each timing offset is replaced with its average over the plurality of times. In some of these embodiments, the average is computed over a predetermined number of times.
• In highly-preferred embodiments of the inventive automatic LAT-determining method, monitoring the timing stability of the first reference channel signal includes monitoring multiple timing offsets LAT₁(X) where X represents the channels with which timing offsets with the first reference channel are computed. Some of these embodiments further include computing a signal characteristic for the plurality of reference channels and determining therefrom which one or more channels among these reference channels has/have not lost timing stability. Some embodiments also include selecting the second reference channel signal from channels which have not lost timing stability, and in some of these embodiments, selecting the second reference channel signal from channels which have not lost timing stability includes computing signal quality.
• In some preferred embodiments, computing a signal characteristic includes computing the frequency content of the signal. In some of these embodiments, computing the frequency content of the signal includes computing a fast Fourier transform (FFT) for a predetermined time period of the signal. In some such embodiments, the computed signal characteristic is the first moment of the signal power determined from the computed fast Fourier transform.
• In some highly-preferred embodiments, computing frequency content of a signal includes segmenting the signal into a plurality of time-overlapping segment signals. In some of these embodiments, weightings are applied to each of the segment signals. In some such embodiments, computing the fast Fourier transform of the signal includes (a) computing a signal-segment fast Fourier transform for each segment signal and (b) averaging each such signal-segment fast Fourier transform to form the fast Fourier transform of the signal. In some of these embodiments, the computed signal characteristic is the first moment of the signal power determined from the fast Fourier transform of the signal.
• In other embodiments, computing the frequency content of the signal includes computing a Haar transform for a predetermined time period of the signal, and in some such embodiments, the computed signal characteristic is the first moment of the signal power determined from the computed Haar transform. In some of these embodiments, the signal is segmented into a plurality of substantially-sequential segment signals. In some embodiments, computing the Haar transform of the signal includes (a) computing Haar transform coefficients for each segment signal, (b) computing absolute values of the coefficients, (c) computing a set of frequency-selective aggregate magnitudes for each segment signal by summing signal-segment Haar transform coefficients having like time scales, and (d) averaging the sets of frequency-selective aggregate magnitudes to form a single set of frequency-selective aggregate magnitudes for the signal. In some of these embodiments, the computed signal characteristic is the first moment of the signal power determined from the frequency-selective aggregate magnitudes.
• In certain other embodiments of the inventive automatic LAT-determining method, the computed signal characteristic is fraction of time within a predetermined time period of the signal at which the absolute value of the signal is above a predetermined threshold.
• In certain other embodiments, the computed signal characteristic is the maximum signal amplitude minus the minimum signal amplitude within a predetermined time period of the signal.
BRIEF DESCRIPTION OF THE DRAWINGS
• FIG. 1 is a schematic block diagram of a method for measuring parameters of MCCE signals, including intracardiac cycle lengths and local activation times and estimates of signal and measurement quality. The steps of the method as illustrated in the block diagram ofFIG. 1 are further detailed in several other schematic block diagrams.
• FIG. 2 is a schematic block diagram illustrating steps of a method to generate absolute-value velocity data of a selected digitized MCCE signal.
• FIG. 3A is an illustration of the operation of the filtering which occurs by applying two differenced sequential boxcar filters to a digitized signal.
• FIG. 3B depicts the absolute value of the output signal of the filtering operation illustrated inFIG. 3A.
• FIG. 4A is a schematic block diagram of the process of determining activations (activity triggers) in the absolute-value velocity signal from an MCCE channel. The steps of this process are applied to more than one channel signal in the method.
• FIG. 4B illustrates the process of identifying activations in an example absolute-value velocity channel signal as processed by the process ofFIG. 4A.
• FIG. 5 is a schematic block diagram of the process of determining the ventricular-channel cycle length, herein called “pulse interval” when associated with a ventricular channel.
• FIGS. 6A and 6B together are a schematic block diagram of the process of determining the reference-channel cycle length in the embodiment ofFIG. 1.
• FIGS. 7A and 7B together are a schematic block diagram of the process of determining local activation time (LAT) for a single mapping point in the embodiment ofFIG. 1.
• FIG. 7C is schematic block diagram of an alternative embodiment to determine LAT for a single mapping point, using additional fiducial times within a reference-channel signal.
• FIG. 8A is a set of MCCE signal plots illustrating an example of the process of determining LAT for a single mapping point as shown inFIGS. 7A and 7B.
• FIG. 8B is a table which illustrates the process by which a specific mapping-channel activation is selected for the determination of LAT for the example ofFIG. 8A.
• FIG. 8C-1 throughFIG. 8C-4 is a set of plots illustrating in detail a selected mapping-channel activation and its corresponding portion of the reference-channel signal which, as illustrated inFIGS. 8A and 8B, are used to determine LAT for a single mapping point.
• FIG. 8D is a table illustrating an embodiment of a method to assess measurement confidence in the method ofFIG. 1, using the examples ofFIG. 8A throughFIG. 8C-4.FIG. 8D also illustrates a second alternative method embodiment to determine an LAT value for a single mapping point.
• FIG. 9 is a schematic diagram illustrating the inclusion of automatic selection of the ventricular and reference channels in the automatic method of measuring parameters of MCCE signals.
• FIG. 10A is a schematic block diagram of the process of automatically selecting a ventricular channel from a set of candidate MCCE channels, specifically illustrating the determination of parameters for a single candidate ventricular channel.
• FIG. 10B is a schematic block diagram of the process of automatically selecting a ventricular channel from a set of candidate MCCE channels, specifically illustrating the automatic selection from among candidate ventricular channels which have had parameters determined in the automatic process ofFIG. 10A.
• FIG. 11A is a schematic block diagram of the process of automatically selecting a reference channel from a set of candidate MCCE channels, specifically illustrating the determination of parameters for a single candidate reference channel.
• FIG. 11B is a schematic block diagram of the process of automatically determining the variability parameter for a single candidate reference channel as used within the process illustrated inFIG. 11A.
• FIG. 11C is a schematic block diagram of the process of automatically selecting a reference channel from a set of candidate MCCE channels, specifically illustrating the automatic selection from among candidate reference channels which have had parameters determined in the automatic process ofFIGS. 11A and 11B.
• FIG. 12 is a matrix which schematically illustrates a series of reference channels and mapping channels among a set of MCCE signals which in an aspect of the inventive method may be processed in parallel to generate multiple LAT maps by various combinations of reference and mapping channels.
• FIG. 13 is a schematic block diagram illustrating an alternative embodiment of the monitoring of cardiac channel quality. The alternative embodiment replaces a portion of the schematic block diagram ofFIG. 1.
• FIG. 14 is a schematic block diagram illustrating an alternative embodiment of the channel selection method illustrated inFIG. 9, adding elements of the automatic channel selection steps within initialization to the real-time operation of the inventive method for measuring parameters of MCCE signals such that cardiac channels may be replaced when the quality of a channel signal degrades during operation of the inventive method.
• FIG. 15 is a high-level schematic block diagram illustrating the steps of an embodiment of the inventive method for transforming LAT values using a second reference channel signal when the timing stability of a first reference channel signal degrades below a stability standard, in order to avoid substantial loss of LAT values in spite of the loss of timing stability.
• FIG. 16 is a table showing exemplary values for timing offsets computed in one step within the embodiment of the inventive method for determining local activation time shown inFIG. 15. The table ofFIG. 16 also presents an example determination of whether or not a loss of timing stability has occurred in the example of the method embodiment ofFIG. 15.
• FIG. 17 is schematic block diagram illustrating the steps of a method for computation of a signal characteristic for use within the method embodiment ofFIG. 15. This signal characteristic computation generates an FFT-based (fast Fourier transform-based) parameter called herein the epoch center-of-power frequency.
• FIG. 17A is an exemplary six-second epoch of a representative cardiac channel electrogram signal.
• FIG. 17B is a plot illustrating one embodiment of weightings used to divide signal epoch data samples into overlapping segments of data. The weightings are applied to the signal data ofFIG. 17A for use within the FFT-based signal characteristic computation ofFIG. 17.
• FIGS. 17C-17G are five plots illustrating the resulting segments having weightings applied to the cardiac channel signal ofFIG. 17A, to be used in the computation of an epoch center-of-power frequency as illustrated inFIG. 17.
• FIGS. 17H-17L are five plots illustrating the segment spectra computed with a fast Fourier transform of the five weighted segment signals ofFIGS. 17C-17G.
• FIG. 17M is a plot of the average signal spectrum of the five segment spectra ofFIGS. 17H-17L and from which an epoch center-of-power frequency is computed.
• FIG. 17N is a table illustrating a method embodiment of the determination of which particular channel(s) have caused the loss of timing stability which was determined to have occurred in the example ofFIG. 16. The signal characteristic used in this example embodiment is the FFT-based signal characteristic computation alternative ofFIG. 17.
• FIG. 18 is schematic block diagram illustrating a first alternative method for computation of a signal characteristic for use within the method embodiment ofFIG. 15. This alternative signal characteristic computation generates a signal characteristic called epoch activity duration.
• FIG. 18A is a plot illustrating the application of an absolute-value velocity filter to the exemplary six-second epoch cardiac channel electrogram signal ofFIG. 17A.
• FIG. 18B is a plot illustrating the computation of an activity duration signal characteristic for the absolute-value velocity epoch signal ofFIG. 18A.
• FIG. 19 is schematic block diagram illustrating a second alternative method for computation of a signal characteristic for use within the method embodiment ofFIG. 15. This alternative signal characteristic computation generates a signal characteristic epoch peak-to-peak amplitude.
• FIG. 19A is a table illustrating the second alternative embodiment ofFIG. 19, showing the peak-to-peak determination of the exemplary six-second epoch of a representative cardiac channel electrogram signal ofFIG. 17A.
• FIG. 20 is schematic block diagram illustrating a third alternative method for computation of a signal characteristic for use within the method embodiment ofFIG. 15. This alternative signal characteristic computation generates a Haar-transform-based parameter called the epoch center-of-power frequency.
• FIGS. 20A-20C are three plots of segments of the exemplary six-second epoch of a representative cardiac channel electrogram signal ofFIG. 17A.
• FIGS. 20D-20F are three plots of the Haar transformation coefficients for the three data segments ofFIGS. 20A-20C, respectively, resulting from the Haar transformation of such three segments of data.
• FIG. 20G is a table detailing the computation of a Haar transformation of a cardiac electrogram signal consisting of 2,048 signal values and resulting in 2,048 Haar transformation coefficients H_i.
• FIG. 20H is a table detailing the computation of a set of eleven frequency-selective aggregate magnitudes A_ifrom the 2,048 Haar transformation coefficients H_iofFIG. 20G.
• FIGS. 20I-20K are three plots of the absolute values of the Haar transformation coefficients shown inFIGS. 20D-20F for the three data segments ofFIGS. 20A-20C, respectively.
• FIGS. 20L-20N are three bar charts of the eleven frequency-selective aggregate magnitudes A_ifor each of the three segments of signal data ofFIGS. 20A-20C.
• FIG. 20P is a bar chart presenting the eleven frequency-selective aggregate magnitudes used in the determination of a COP-frequency signal characteristic using the alternative method ofFIG. 20.
DETAILED DESCRIPTION OF PREFERRED EMBODIMENTS
• FIG. 1 illustrates one embodiment of a method for measuring parameters of multi-channel ECG signals.FIG. 1 is a high-level schematic block diagram of a method which measures intracardiac cycle lengths and local activation times on a near-real-time basis as part of a system to generate maps (e.g., computer displayed 3D presentations of the distribution of voltages and activation times across cardiac structures) and provides feedback regarding signal quality and measurement confidence.
• Several other figures in this document relate to the method ofFIG. 1, and the steps presented in the schematic block diagrams of these other figures are nested within the high-level schematic block diagram ofFIG. 1, as will be described below. In addition, the method also includes initial channel selection steps which occur prior to the steps ofFIG. 1. These are illustrated and described later in this document, inFIGS. 9 through 11C.
• Referring toFIG. 1, anembodiment10 of the method includes a flow loop of method steps which is initiated by arequest12 to map a point, and each time a mapping-point request12 is generated, the method proceeds through the steps shown inFIG. 1. The flow chart element labeled withreference number14 indicates that the flow loop waits to receiverequest12. During a procedure in which the method is used, an electrophysiologist (EP doctor) is maneuvering an electrode-tipped catheter (mapping catheter) through and around the chambers, arteries and veins of a patient's heart. The electrode on this maneuvered catheter provides the mapping-channel signal. When the EP doctor determines that the maneuvered catheter electrode is in a desired position, the EP doctor activates a signal asrequest12 to map a point. A plurality of map points constitute the map.
• Generating the map during this procedure involves time measurements made between the MCCE signals of the mapping electrode and a reference electrode. (As used herein, electrodes are positioned to provide signals to channels. Thus, for example, the mapping electrode provides the signal for the mapping channel.) The reference electrode is positioned before mapping begins in a location that is expected to remain constant during the mapping process and that will generate stable and repetitive electrical signals.
• Each electrode develops an electrical signal when muscle cells in contact with the electrode change their cell membrane potentials. These electric potentials change as the cells mechanically contract. Nerve cells, which do not contract, also can be in contact with electrodes and produce electrical signals.
• The map being generated represents a particular heart rhythm being studied, such as tachycardia. The reference-channel and mapping-channel signals are both cyclical and have substantially the same cycle length (CL). The reference-channel signal represents a zero-phase or index moment of the particular cardiac cycle, and the local activation time (LAT) measurements (time difference between mapping and reference-channel signals) indicate the sequence of muscle and nerve cell activation of various points (map points) in the cardiac structure. This time sequence and its physical course around the anatomy of the heart are the information the EP doctor needs to determine how to apply therapy. The term “local” refers to the fact that the measurement applies to the heart cells in contact with the electrode and to signals with respect to a reference-channel signal, and this information is translated to a position on a three-dimensional (3D) image of the heart chamber.
• Activation time is measured relative to one or more activations at the reference electrode and may be positive or negative. A local activation time which is negative by more than a half of one cycle length may also be recognized as being positive at a corresponding time less than a half of one cycle length. Local activation times may be defined as being relative to the nearest activation in the reference channel.
• Positioning of the mapping catheter is guided at times by fluoroscopic imaging. At a position of interest, the EP doctor generatesrequest12 to trigger the system to make measurements from the MCCE signals available from the maneuvered catheter and other more stationary catheters and body surface electrodes. These measurements at mapping points are represented graphically, usually by color, on a 3D image of the heart chamber of interest. These points may be requested at irregular intervals of several seconds to perhaps minutes, depending on when the EP doctor maneuvers the mapping catheter to a point at which measurements should be taken.
• Whenrequest12 is received, measurements are made using an “epoch” of the most recent 6 seconds of MCCE signals. Inembodiment10, the 6-second length of this epoch should not be taken as limiting. The epoch is a preset time window of MCCE signals, and its 6-second length is chosen here inembodiment10 such that selected signals during the preset time window contain a suitable number of electrical events to permit the analysis to be performed. During such mapping procedure, at least one mapping channel and at least one reference channel are used. At some points withinembodiment10, as will be described later in this document, the epoch is divided into three equal periods of time, and six seconds is chosen here since a 2-second period will almost always contain at least one heart beat (or cell activation) for all heart rates above 30 beats per minute.
• As the mapping catheter is moved, it is important that its electrode be in place at the selected location for a period of time (dwell time) long enough to obtain a suitable signal. Inembodiment10, such dwell time is about 2 seconds. Thus, whenrequest12 is received, the epoch consists of 6 seconds of data on other channels being used and 2 seconds of data on the mapping channel. (The 6 seconds of data may consist of the immediate past 4 seconds of the data plus 2 seconds of data generated afterrequest12 occurs. The 6 seconds of data in an epoch may also be the 6 seconds of data immediately preceding therequest12, since it may be that the mapping catheter has already been in a stable position for the 2 seconds prior to the triggering ofrequest12. Other possible strategies for acquiring the epochs of data are also possible.)
• In the high-level schematic block diagram ofFIG. 1, afterrequest12 is received, ending the wait inmethod step14, a determination of the pulse interval in the ventricular channel is performed inmethod step15. Details of ventricular pulse-interval determination15 are detailed in the schematic block diagrams and example signals ofFIG. 2 throughFIG. 5, all of which will be described later in this document.
• Following ventricular pulse-interval determination15, adetermination16 of the intracardiac cycle length in the reference channel is performed. (Method step16 is shown inFIG. 1 as determining tachycardia cycle length since the method is intended primarily for monitoring cardiac parameters in the treatment of patients in tachycardia. The use of the term “tachycardia” is not intended to be limiting. The method is applicable to measurement of all types of cardiac arrhythmias as well as normal heart rhythms.) Details of intracardiaccycle length determination16 are detailed in the schematic block diagrams and example signals ofFIGS. 2 through 4B andFIG. 6, all of which will be described later in this document. Intracardiaccycle length determination16 depends on the value of ventricular pulse interval established indetermination15.
• Decision step18 followsdetermination16 such that the cycle length determined instep16 is compared to a cycle-length-change criterion indecision step18, and if the cycle length has not exceeded the cycle-length-change criterion, the method proceeds. If, however, the cycle-length-change criterion is exceeded, the EP doctor is alerted inmethod step20 in order that steps may be taken by the EP doctor during the mapping procedure to evaluate the impact of such a change.
• A cycle-length-change criterion applied inmethod step18 may be based on an absolute time difference in cycle length from a previous cycle length or on the average of a plurality of previous cycle lengths. Or it may be based on a percentage change from such quantities. One useful previous cycle length is the initial or starting cycle length of the reference channel, established at the beginning of the mapping procedure. A local activation time map is related to a particular rhythm so that if there is too great a change in cycle length, the EP doctor may choose to start a new map, or in fact may determine that mapping is no longer appropriate at such time. A value for the percentage change which triggers an alert inmethod step20 may be that the current reference-channel cycle length (determined in method step16) is found to differ from the starting cycle length by more than 10%. Such value is not intended to be limiting; other values may be found to provide adequate warning to the EP doctor.
• Embodiment10 of the method then proceeds to acomputation22 of the local activation time (LAT) associated with the map point being analyzed. Details of localactivation time computation22 are detailed in the schematic block diagram ofFIGS. 7A-7C which will be described later in this document, and examples of such determination are illustrated inFIGS. 8A through 8D.
• Embodiment10 of the method for measuring parameters of MCCE signals includes steps forevaluation24 of signal quality andevaluation26 of measurement confidence, both of which are applied withinembodiment10 to monitor the measurement process. In each case, that is, reduced signal quality as determined instep24 and reduced measurement confidence instep26, the EP doctor is alerted (user alerts28 and30, respectively) that such conditions have been detected. One embodiment of a method to measure signal quality inmethod step24 is included in the steps illustrated inFIG. 4A and will be discussed later in this document. One embodiment of a method to assess measurement confidence inmethod step26 is illustrated in the example ofFIG. 8D described later in this document.
• As shown inFIG. 1, the method ofembodiment10 provides (in step32) the map point and its related measurement data to a computer at least for display to the EP doctor during the procedure and for storage in memory for later analysis. The system then returns vialoop path33 to wait for the nextmapping point request12 atstep14.
• FIGS. 2, 3A and 3B illustrate an embodiment of a portion of the steps of the method further detailed inFIGS. 4A-11C.FIG. 2 is a schematic block diagram illustrating detailed steps withinembodiment10 by which a selected digitizedMCCE channel signal36iis filtered to generate a corresponding absolute-value velocity signal36o. The steps ofFIG. 2 are applied to various signals withinmethod embodiment10, as indicated later in the description below.
• InFIG. 2, the combined steps, low-pass filter38, first-difference filter40, and absolute-value filter42, are together shown as an absolute-value velocity filter34. The first two steps of absolute-value velocity filter34, low-pass filter38 and first-difference filter40, are together shown as abandpass filter44 which generates a filteredvelocity signal41 of an input signal.
• As shown inFIG. 2, aninput signal36iis a 6-second preset time window (epoch) of digitized data from a selected MCCE signal. Low-pass filter38 operates oninput signal36ifollowed by first-difference filter40, and together these two filters generate a digital stream ofdata41 which corresponds to the filtered velocity (first derivative) ofinput signal36iwith certain low and high frequencies filtered out. That is,filter44 is a bandpass filter. Absolute-value filter42 simply applies an absolute-value operation (rectification) to filteredvelocity signal41 from first-difference filter40 to generate output signal36owhich is an absolute-value velocity signal ofinput signal36i.
• One embodiment of applying acombination44 of low-pass filter38 and first-difference filter40 to a digitized signal is what is called herein “two differenced sequential boxcar filters,” and such filtering embodiment is illustrated inFIG. 3A, in which an example digitized signal46 is shown both graphically (46g) and numerically (46n). Seven pairs of “boxcars”48billustrate the sequential operation ofboxcar filter48.
• Referring toFIG. 3A, each pair ofboxcars48binboxcar filter48 is four time samples in length, and twoboxcars48bare such that one follows the other immediately in time. (Only two of the 14boxcars48bare labeled.) The sum of the four time samples of digitized signal46 in eachboxcar48bis calculated. Thus, for example, theleft boxcar48bof the uppermost (first in time) pair as shown holds the sum of the four time samples it subtends, and theright boxcar48bof this pair holds the sum of the four time samples it subtends. These two sums are 7 and 12, respectively, and the difference between the right boxcar value and the left boxcar value is 12−7=5. Thisdifferenced value 5 is shown to the right of theuppermost boxcar48bpair, and seven such values, indicated byreference number50, are shown to the right of the seven examplesequential boxcar48bpairs. Thisoutput signal50 is shown both numerically as50nand graphically as50g.Filter output signal50 is shown for the seven time samples between the dotted lines labeled52aand52b.
• In the example ofFIG. 3A, eachboxcar48bhas a boxcar-width w_Bof four samples. The value of w_Bdetermines the frequency response ofboxcar filter48, or the amount of smoothing provided byboxcar filter48. Larger values of w_Bproduce a lower central frequency ofboxcar filter48 and therefore more smoothing of the signal on which it operates. Such relationships are well-known to people skilled in the art of digital filtering. Any specific value for w_Bused herein is not intended to be limiting. However, for the embodiments exemplified herein, it has been found that values of w_Bof around 4 are appropriate for use on intracardiac signals, and values of w_Baround 20 are appropriate for ventricular channels. For MCCE signals digitized every 1 millisecond and for sequential four-sample-long boxcar filters48 (w_B=4) illustrated inFIGS. 3A, 3B and 5, the resulting bandpass filter has a center frequency of 125 Hz.
• The operation of the two differenced sequential boxcar filters48 performs low-pass filtering and differentiation to input signal46 such thatfilter output50 is proportional to the velocity of bandpass-filtered digitized signal46. No scaling has been applied in this example, but such lack of scaling is not intended to limit the meaning of the term two differenced sequential boxcar filters.
• FIG. 3B, shown below and to the left ofFIG. 3A, simply graphically illustrates the absolute value ofoutput signal50 as exampled inFIG. 3A. The absolute value ofoutput signal50 is processed by absolute-value filter42 as shown inFIGS. 3A and 3B. This absolute-value velocity signal is output signal36oofFIG. 2.
• Some steps of the method as illustrated inembodiment10 include the identification of activations or activity triggers within one or more channel signals of MCCE signals. Activations (activity triggers) are the electrical activity associated with the initiation of the depolarization of the heart muscle cells which occurs during a heart beat, progressing like a wave through the various portions of the cardiac structure and causing the heart to pump.
• FIG. 4A is a schematic block diagram of aprocess58 of determining activations (activity triggers) in an absolute-value velocity signal. The steps ofprocess58 may be applied to more than one channel signal.
• In the embodiment ofFIG. 4A, asignal60 which is 6 seconds in duration (6-sec epoch) and is the absolute-value velocity of an MCCE channel signal, is divided into three 2-second “chunks” inmethod step62. In method steps64-1,64-2 and64-3, these three chunks are processed to find three signal maxima (max1, max2, and max3), one for each of the three signal chunks. These three values (max1, max2, and max3) are inputs tomethod step66 which selects the maximum MAX among the three inputs andmethod step68 which selects the minimum MIN among the three inputs. The values MAX and MIN are in turn inputs tomethod step70 which determines an estimate SI for signal irregularity. Inmethod step70, signal irregularity SI is estimated as SI=MAX−MIN. A larger difference between the maximum (MAX) and minimum (MIN) values of the chunk maxima (max1, max2, and max3) indicates that there is more irregularity among the heart beats withinepoch60 being processed. Signal irregularity SI is related to the variations in the “shape” of the activations in MCCE signals while other measurements described later in this document relate to variations in the time of activations.
• The value MIN represents an estimate SS of signal strength. SS is multiplied by 0.5 (threshold factor) inmethod step72 to determine a value for an activation threshold AT to be used instep74 to determine the occurrence of activations within the MCCE signal being processed. The value (0.5) of the threshold factor applied inmethod step72 of this embodiment is not intended to be limiting. Other values for the threshold factor maybe be applied in embodiments of the method.
• Signal irregularity SI and signal strength SS are used in conjunction with an estimate of signal noise N_Sto provide an estimate of signal quality SQ inmethod step79. Inmethod step78, signal60 (provided byflow path60a) is processed to compute its median over the entire six-second epoch, and such median is multiplied by 2 to produce estimate N_Sof signal noise. Inmethod step78, the calculation of the median ofsignal60 may be done using a normal median or a set-member median. For such large data sets (e.g.; 6 seconds at 1,000 samples per second), it has been found that using the set-member median is computationally convenient and highly suitable. Instep79, signal quality SQ is computed as SQ=SS−SI−2N_S.
• The factor of 2 applied inmethod step78 and the factor of 2 applied inmethod step79 are both not intended to be limiting. Other values for such factors may be used. The size of the factor instep78 is related to ensuring that the estimate of noise N_Sinsignal60 is a good representation of the noise level insignal60. The size of the factor instep79 is related to the relative weight given to noise estimate N_Scompared to those given to signal strength SS and signal irregularity SI in generating the estimate for signal quality SQ. The values of 2 for both of these factors have been found to provide good performance for estimating noise N_Sand signal quality SQ.
• FIG. 4B illustratesmethod step74 ofFIG. 4A, the process of identifying activations in an example absolute-valuevelocity channel signal60 as processed by the method steps ofFIG. 4A. The signal ofepoch60 being processed is an input tomethod step74 as indicated by thesignal flow path60a. A portion ofexample epoch60 is illustrated inFIG. 4B. Activation threshold AT is shown as a dotted line AT parallel to the time axis and intersectingsignal60 at points76. (Eleven signal crossings are shown; one such point is labeled76a, one is labeled76b, and one is labeled76c).
• As indicated inmethod step74 ofFIG. 4A, activations inepoch60 being processed are indicated by identifyingthreshold crossings76 before which signal60 does not cross activation threshold AT for at least T_BTmilliseconds. The value of before-threshold time T_BTchosen may vary according to the type ofMCCE signal60 being processed. For example, it has been found that T_BT=120 msec is an appropriate value when a ventricular channel is being analyzed, and that T_BT=90 msec is an appropriate value when an intracardiac channel is being analyzed. These values for T_BTare not intended to be limiting; the selection of a value for T_BTis based on choosing a value by which a reliable differentiation between subsequent activations and amongthreshold crossings76 within an individual activation can be achieved.
• In the example ofFIG. 4B, the activation labeled75 shown includes sixthreshold crossings76 as indicated by dotted circles occurring in rapid succession, the first being threshold crossing76band the last being threshold crossing76c. A portion of aprevious activation77 withinsignal60 is also shown inFIG. 4B. Inactivation77, fivethreshold crossings76 occur in rapid succession, the last of which is labeled76a.
• The time difference between threshold crossing76aassociated withactivation77 and threshold crossing76bassociated withactivation75 is about 185 msec as shown inFIG. 4B. In this example, 185 msec is longer than either of the example values for T_BT; thus threshold crossing76bis determined to be the leading edge ofactivation75, and the time at which threshold crossing76boccurs is determined to be activation time t_ACT. In this example,threshold76bis the only such threshold crossing illustrated inFIG. 4B.
• FIG. 5 is a schematic block diagram of anembodiment80 of a process of determining the ventricular-channel pulse interval (heart rate). In the automatic method of measuring parameters of MCCE signals, ventricular-channel pulse-interval information is used in the calculations to determine intracardiac-channel parameter calculations. The steps ofembodiment80 ofFIG. 5 analyze an absolute-value velocity ventricular-channel signal epoch82, again 6 seconds in duration. Inmethod step84, activations withinepoch82 are identified by applyingsteps58 as illustrated inFIG. 4A. As indicated in method steps82 and84, when processing a ventricular channel, values for boxcar-width w_Band before-threshold time T_BTmay differ from those used for intracardiac channels. In the embodiment ofFIG. 5, these values are w_B=20 and T_BT=120 msec.
• Activations identified inmethod step84 each have an activation time and for purposes of description, there are n such activation times. Inmethod step86, all activation intervals I_iare computed. There are n−1 activation intervals I_icomputed as follows:
• $\hspace{1em}\begin{array}{c}{I}_1=t_2-{\mathrm{<figure-callout\; label="t"\; filenames="US20160089048A1-20160331-D00003.png,US20160089048A1-20160331-D00007.png"\; state="\{\{state\}\}">t</figure-callout>}}_1\\ •\\ I_i=t_{i+1}-{\mathrm{<figure-callout\; label="t"\; filenames="US20160089048A1-20160331-D00003.png,US20160089048A1-20160331-D00006.png"\; state="\{\{state\}\}">t</figure-callout>}}_2\\ •\\ I_{n-1}=t_n-t_{n-1}\end{array}$
• Inmethod step88, a maximum interval MAX_PIof the n−1 activation intervals I_iis computed, and instep90, the minimum interval MIN_PIof the n−1 activation intervals I_iis computed. Inmethod step92, a range R_PIfor activation intervals I_iis computed as the difference between MAX_PIand MIN_PI.
• The n activation times t_iare also used inmethod step94 to compute all double-intervals D_iof ventricular-channel signal epoch82. There are n−2 double-intervals D_i, and such double-intervals D_iare computed as follows:
• $\hspace{1em}\begin{array}{c}{\mathrm{<figure-callout\; label="D"\; filenames="US20160089048A1-20160331-D00003.png,US20160089048A1-20160331-D00007.png"\; state="\{\{state\}\}">D</figure-callout>}}_1=t_3-{\mathrm{<figure-callout\; label="t"\; filenames="US20160089048A1-20160331-D00003.png,US20160089048A1-20160331-D00007.png"\; state="\{\{state\}\}">t</figure-callout>}}_1\\ •\\ D_i=t_{i+2}-{\mathrm{<figure-callout\; label="t"\; filenames="US20160089048A1-20160331-D00003.png,US20160089048A1-20160331-D00006.png"\; state="\{\{state\}\}">t</figure-callout>}}_2\\ •\\ D_{n-2}=t_n-t_{n-2}\end{array}$
• Inmethod step96, the normal median M_DIof all double-intervals D_iis computed, and instep98, the estimate PI of ventricular-channel pulse interval is computed as
• 
  PI=M_DI/2
• Thus, method steps ofprocess80 generate an estimate of ventricular pulse interval PI and provide an estimate of the range R_PIover which ventricular pulse interval PI varies. The value of pulse interval PI is used in the determination of reference-channel and mapping-channel cycle lengths and is reported as a heart rate HR for the patient being monitored. Heart rate HR in beats per minute (bpm) is determined inmethod step99 from pulse interval PI (in msec). (For computational convenience instep96, a set-member median calculation may be used in place of the normal median calculation.)
• FIGS. 6A and 6B together are a schematic block diagram of aprocess220 of determining the reference-channel cycle length CL. (Method step16 ofFIG. 1 includesprocess220 along with other elements of the reference-channel cycle-length determination as illustrated in the embodiments ofFIGS. 2 through 4B.) The steps ofFIGS. 6A and 6B process an absolute-value velocity reference-channel signal epoch222, again 6 seconds in duration. Instep224, reference-channel signal strength SS, signal irregularity SI, and noise N_Sare computed using the steps ofFIG. 4A. Method steps72ands74 ofFIG. 4A are not used in the processing of reference-channel absolute-value velocity signal222; individual activations withinsignal222 are not detected. Signal strength SS, signal irregularity SI, and noise N_Sare used to determine signal quality SQ, and noise N_Sis used in the computation of the magnitude-coincidence autocorrelation function ACF inmethod step226.
• Inmethod step226, a magnitude-coincidence autocorrelation is performed on the data in absolute-value velocity reference-channel signal epoch222. (The computed autocorrelation function is indicated by the term ACF.) The threshold value for the magnitude-coincidence autocorrelation is dependent on noise N_Sinsignal222 as described in the summary section above which defines magnitude-coincidence autocorrelation. As applied inmethod step226, the value of the threshold T_ACis set to ensure that the thresholding process selects events which are significant events withininput signal222. In one embodiment,
• 
  T_AC=2·N_Swhere noiseN_S=2·(median(input)+1).
• The “1” is added to the median for computational convenience and to avoid singular conditions within the system. Values other than 1 may be used and other ways to set threshold T_ACmay be used; this specific expression for T_ACis not intended to limit the scope of this method.
• The remaining method steps ofprocess embodiment220 inFIGS. 6A and 6B, steps228-252, encompass analysis of ACF to determine reference-channel cycle length CL. As described above, ACF is a function of lag, and in this analysis, values of lag in ACF are analyzed with respect to ventricular-channel pulse interval PI such that reference-channel cycle-length CL is determined inprocess220 based on estimates of pulse interval PI.
• Inmethod step228, a minimum of ACF at values of lag less than about 200 msec is identified. (200 msec is a preset lag threshold.) The lag at this minimum in ACF is labeled W and is an estimate of activity width. The lag threshold value of 200 msec for searching for activity width is chosen such that the width of activations expected for most intracardiac-channel signals will be found at lag values less than 200 msec. The search window (preset lag threshold) should be shorter than the shortest expected value of reference-channel cycle length and longer than the width of activations in the reference-channel signal. Since activations typically are significantly shorter than CL, it is straightforward to set the range to an appropriate value. 200 msec has been found to be a useful value. However, the specific value of 200 msec for the preset lag threshold is not intended to be limiting.
• Inmethod step230, the maximum peak P₁is found in ACF for values of lag greater than W; CL is set at the value of lag CL₁where ACF has its maximum peak P₁for lag greater than W; and an interim peak amplitude P_CLis set to P1. (P_CL, P₁, P₂, P₃, CL₁, CL₂and CL₃are interim values in the steps ofprocess220.) Then inmethod step232, if CL₁is very near (within ±20 msec) double the ventricular pulse interval PI, then process220 proceeds tomethod step234. If CL₁is not very near 2PI, then process220 proceeds tomethod step242 inFIG. 6B. (Points A and B inprocess220 represent common points joiningFIGS. 6A and 6B.) The situation of CL₁being very near double the ventricular pulse interval PI may occur in a condition in which a slight alternation of ventricular pulse interval PI in a pattern of bigeminy (long, short, long, short, etc.) causes the autocorrelation to be slightly stronger at the repeating double interval. Consecutive single intervals are only slightly different, but magnitude-coincidence autocorrelation is very sensitive to the slight timing differences.
• Throughoutprocess embodiment220 of determining reference-channel cycle length CL, there are several time intervals which are used to identify certain values in ACF such as the ±20 msec “nearness” criterion inmethod step232. These occur in method steps232,234,238,242,248, and250. In each such occurrence, these specific values have been found to perform well in the embodiment ofprocess220. (The “nearness” criteria are also referred to as lag intervals. The lag intervals in the method steps ofprocess220 are:step232, a first lag interval;step234, a second lag interval;step238, a third lag interval;step242, a fourth lag interval; andsteps248 and250, a fifth lag interval.)
• Inmethod step234, the maximum amplitude P₂of ACF is identified within a lag interval of ±40 msec of ventricular pulse interval PI, CL₂is set to the value of lag at maximum P₂, andprocess220 proceeds tomethod step236. Inmethod step236, if the amplitude P₂is greater than half of peak amplitude P_CLand if, inmethod step238, CL₂is within 20 msec of CL₁/2, then inmethod step240, CL is set to CL₂, P_CLis set to P₂, andprocess220 proceeds to step242. If both of these two conditions (insteps236 and238) are not true,process220 proceeds to step242 without setting CL to CL₂and P_CLto P₂. Method step238 distinguishes peak P₂from a maximum on one of the boundaries of the ±40 msec lag interval inmethod step234. If the P₂is not greater than half of P₁, then the process proceeds tomethod step242.
• Inmethod step238, if CL₁/2 is within 20 msec of CL₂, then CL is set to CL₂inmethod step240 and the process proceeds tomethod step242. If CL₁/2 is not within 20 msec of CL₂, then the process proceeds tomethod step242 without setting CL to CL₂.
• Inmethod step242, if CL (set inmethod step230 or method step240) is within 60 msec of ventricular pulse interval PI, then process220 proceeds tomethod step244. If CL is not within 60 msec of PI, then the process ends and the reference-channel cycle length is either CL=CL₁as set inmethod step230 or CL₂as set inmethod step240.
• Inmethod step244, the maximum amplitude P₃of ACF is identified within the lag interval between lag=CL/6 and lag=2CL/3, interim value CL₃is set to the lag at maximum P₃, andprocess220 proceeds tomethod step246. Inmethod step246, if amplitude P₃is greater than half of the amplitude P_CLat CL (CL is either the lag CL₁at peak P₁or the lag CL₂at peak P₂), then process220 proceeds to method steps248 and250. If the amplitude P₃does not satisfy the criterion inmethod step246, then process220 ends and the value of reference-channel cycle length CL is either CL=CL₁as set inmethod step230 or CL₂as set inmethod step240.
• It is possible that there may be a significant peak in ACF between lag=CL/6 and lag=2CL/3. Method steps248 and250 are parallel steps which, if either of the criteria in these steps is satisfied,process220 proceeds tomethod step252 in which the reference-channel cycle length CL is set to CL=CL₃andprocess220 ends. If neither of these two criteria is satisfied,process220 ends and reference-channel cycle length CL is either CL=CL₁as set inmethod step230 or CL₂as set inmethod step240. The criteria in method steps248 and250 check whether peak P₃has a value of lag wherein CL is within 20 msec of either 2CL₃or 3CL₃. If either condition is true, then, as stated above, reference-channel cycle length CL is set to CL₃andprocess220 ends. The situation of a proper reference-channel cycle length CL being at ⅓ or ½ of ventricular pulse interval PI is related to 3:1 or 2:1 atrio-ventricular conduction with the artificial enhancement of the ACF peak at pulse interval PI because of the ventricular artifact that occurs for some of the atrial activations.
• The methods just described can be summarized as three distinct and separable strategies. First is the use of the autocorrelation function to identify repeating cycles in the cardiac rhythm with maximum use of all the data available and little dependance on shape, no dependance on threshold-crossing jitter, and robust-to-occasional noise glitches. The second important strategy is avoiding the choice of a false cycle length at twice the ventricular pulse interval because the ventricular response slightly alternates in a pattern of bigeminal timing. The third important strategy is avoiding the choice of a cycle length equal to the ventricular pulse interval because ventricular far-field distortions may occur in atrial signals during 2:1 or 3:1 atrio-ventricular conduction. These three strategies are useful separately but more so in combination.
• FIGS. 7A and 7B together are a schematic block diagram of anembodiment100 of the process of determining local activation time (LAT) for a single mapping point in the method of measuring parameters of MCCE signals.FIG. 7A illustrates three MCCE signals (6-sec epochs) on which computations are performed, as has been described above, in order to provide results which are used in the determination of LAT for a single mapping point. A ventricular-channel epoch102 and a reference-channel epoch108 are coincident in time, and a mapping-channel 2-sec epoch114 is coincident with the last 2 seconds ofepochs102 and108.FIG. 7A includes a legend which defines the terminology used inFIGS. 7A and 7B.
• Inmethod step104, ventricular-channel epoch102 is processed with the steps ofFIG. 4A and produces a set of ventricular-channel activation times t_V-ACTand estimates of signal quality SQ and signal irregularity SI forepoch102. The ventricular-channel activation times t_V-ACTare used in the determination of LAT as indicated by the circle labeled V which is common with the same such circle inFIG. 7B. Inmethod step106, ventricular-channel pulse interval PI and range of pulse intervals R_PIare computed for presentation to the EP doctor using the steps shown inFIG. 5.
• In a similar fashion, inmethod step110, reference-channel epoch108 is processed with the steps ofFIG. 4A and produces estimates of signal quality SQ and signal irregularity SI forepoch108. Within the method steps ofFIG. 4A, a velocity signal for the reference channel is computed, and it is used in the determination of LAT as indicated by the circle labeled R_velwhich is common with the same such circle inFIG. 7B. Inmethod step112, reference-channel cycle length CL is determined using the steps shown inFIGS. 6A and 6B, and reference-channel cycle length CL is used in the determination of LAT as indicated by the circle labeled CL which is common with the same such circle inFIG. 7B.
• Inmethod step114, mapping-channel epoch114 is processed with the steps ofFIG. 4A and produces a set of mapping-channel activation times t_M-ACTand estimates of signal quality SQ forepoch114. (Epoch114 is not sufficiently long to determine a useful estimate of signal irregularity SI. However, if a longer epoch length is used, SI may be estimated inmethod step116.) Mapping-channel activation times t_M-ACTare used in the determination of LAT as indicated by the circle labeled M which is common with the same such circle inFIG. 7B. Within the steps ofFIG. 4A, a velocity signal for the mapping channel is computed, and it is used in the determination of LAT as indicated by the circle labeled M_velwhich is common with the same such circle inFIG. 7B.
• FIG. 7B shows a continuation of the flow chart ofembodiment100. Inputs to the method steps ofFIG. 7B have been computed in the method steps ofFIG. 7A, and these inputs are illustrated by the circles labeled as described above. Inmethod step118, a mapping-channel activation for LAT determination is selected from among the activations and corresponding mapping-channel activation times t_M-ACTdetermined inmethod step116. The selection of such activation instep118 includes the maximization of an activation selection score A_SC, a value for which is computed for each candidate activation among the set of mapping-channel activations. Details ofmethod step118 are described later in this document in the discussion of the example ofFIGS. 8A and 8B.
• After selecting the specific mapping-channel activation to be used to determine LAT inmethod step118, a mapping-channel fiducial time t_Mis found inmethod step120. In determining LAT, a more precise representation of event times is required than the threshold-crossing determination of activation detection inmethod step74. In this document, “fiducial time” is the term used to indicate such a more precise determination of an event (activation) time. “Fiducial time” as used herein represents the instant within an MCCE signal at which a depolarization wavefront passes below the positive recording electrode in either a bipolar or unipolar MCCE signal.
• As is well-known to those skilled in the field of electrophysiology, one good representation of fiducial time is the instant at which a signal exhibits its maximum negative velocity. Thus, one embodiment ofmethod step120 includes determining mapping-channel fiducial time t_Mas the time at which the maximum negative velocity occurs within the selected activation of the mapping channel. In a similar fashion, a reference-channel fiducial time t_Ris found inmethod step122. Reference-channel fiducial time t_Ris the time at which the maximum negative velocity occurs within ±CL/2 of mapping-channel fiducial time t_M.
• The use of the time of maximum negative velocity as the fiducial time is not intended to be limiting. Other indications of precise depolarization event times may be used in determining the fiducial times.
• Inmethod step124, the local activation time LAT for a position at which the mapping-channel electrode is located within the heart is computed as LAT=t_M−t_R. Local activation time LAT is determined relative to the selected reference channel, and values of LAT at a plurality of locations within the region of the heart being mapped are determined during the process of building an LAT map. If the quality of the channel signals being processed degrades before mapping is complete such that mapping cannot be continued, a new map must be generated. Local activation times may be positive or negative times (occurring after or before the corresponding activation event in the reference channel).
• FIG. 7C is a schematic block diagram of analternative embodiment122′ of the process by which an LAT value is determined for a single mapping point. (Embodiment122′ ofFIG. 7C is an alternative embodiment to method steps122 and124 ofFIG. 7B.)FIG. 7C will be described later in this document, after the example ofFIGS. 8A-8D is described.
• FIG. 8A throughFIG. 8D together illustrate in more detail the process of determining LAT for a single mapping-channel electrode location.FIG. 8A is a set of exemplary MCCE signal plots. At the top ofFIG. 8A is a six-second epoch of an ECG reference-channel signal108. At the bottom ofFIG. 8A is a six-second epoch of an ECG ventricular-channel signal102 time-coincident with reference-channel signal108. In the middle and to the right ofFIG. 8A is a 2-second epoch of an MCCE mapping-channel signal114 time-coincident with the final 2 seconds of reference-channel signal108 and ventricular-channel signal102. (Note that inFIG. 8A, the signal traces illustrated are the MCCE voltage signals, not absolute-value velocity signals which are created during signal processing represented in thesteps100 ofFIG. 7A.)
• Also illustrated inFIG. 8A is magnitude-coincidence autocorrelation ACF126 of the absolute-value velocity of reference-channel signal108.ACF126 is used to determine reference-channel cycle length CL in a process such asembodiment220 inFIGS. 6A and 6B.ACF126 inFIG. 8A is annotated to show activity width W. Peak P₁occurs at a lag value of 342 msec, and the lag at peak P₁is the reference-channel cycle length CL in this example.
• Ventricular-channel activations identified inmethod step104 are shown inFIG. 8A. Nine such activations, represented by ventricular-channel activation times128, were identified in ventricular-channel signal102. In a similar fashion, mapping-channel activation times130 are shown inFIG. 8A; four mapping-channel activations were identified inmethod step116. Note that in this example, the first and last activations in the 2-second epoch of mapping-channel signal114 were not identified by the threshold-crossing activation detection process ofmethod step74 ofFIG. 4A. (For description purposes, all six activations in mapping-channel signal114 are labeled132athrough132finFIG. 8A even though only four such activation were detected. The number132 is not repeated with the letters a-f for simplicity in the figure.Activations132aand132fwere not detected.)FIG. 8B presents a table detailing method steps118 by which a specific mapping-channel activation is selected for the determination of LAT for the example ofFIG. 8A.FIG. 8B includes a legend further defining the terms utilized in construction of and computations within the table for this example. Inmethod step112 ofFIG. 7A, reference-channel cycle length CL was determined to be 342 milliseconds, shown inFIG. 8B at the top of the table.
• As mentioned above, local activation time (LAT) is measured by the time difference between a fiducial time t_Min an activation in the mapping channel and its corresponding fiducial time t_Rin the reference channel. As part of this determination, an activation within the mapping-channel signal114 must be selected for such computation, inmethod step118. This selection process includes: (a) for each mapping-channel activation i, determining the time t_NV(i) to the nearest ventricular-channel activation for each mapping-channel activation; (b) for each mapping-channel activation i, determining the deviation D_P(i) from CL of the time to the previous mapping-channel activation i−1; and (c) for each mapping-channel activation i, determining the deviation D_F(i) from CL of the time to the next (future) mapping-channel activation i+1. The mathematical representations of these determinations are shown in the legend ofFIG. 8B.
• To generate a full map of local activation times, often a large number of individual points must be determined. This can be a time-consuming process. It is therefore desirable to determine each individual value of LAT as quickly as possible once a new position of the mapping-channel electrode being manipulated by the EP doctor is established. It has been found that about 2 seconds is often required to make a good determination. At typical intracardiac heart rates being measured, only a few activations occur in the mapping channel during a 2-second epoch period, so it is helpful to increase the number of candidate activations by adapting to situations where an activation is “missing” due to a failed activation detection or to a simple epoch-end timing situation. The method includes a beginning-of-data rule and an end-of-data rule to increase the number of candidate mapping-channel activations. These special rules are as follows:
• Beginning-of-data rule: In some cases, the first detected activity may be very near the beginning of available data. If the expected previous activity to a detected activity would be located before the beginning of the mapping-channel epoch, then there is no evidence that detections failed and the value for D_P(i) for such a candidate activation is presumed to be 0. However, if the amount of time in the available data in the mapping-channel epoch is longer than the expected cycle length CL, then it is likely that an activation failed to be detected due to some kind of noise in the mapping-channel signal, an irregular signal, or an insufficiency in the detection algorithm. In this case, D_P(i) is set to t_M-ACT(i)−CL, but not less than 0, where CL is the reference-channel cycle length.
• End-of-data rule: This rule is symmetrical to the beginning-of-data rule and is created to handle the same available data constraint at the end of the data. WO for only the last candidate mapping-channel activation is set to 0 if the last detected activity is within one reference-channel cycle length CL of the end of data. However, there may be more time in the available mapping-channel epoch data than one CL after the last detected activation. In this case, it is very likely that some kind of noise in the mapping-channel signal, an irregular signal, or an insufficiency in the detection algorithm caused a failed activation detection. In this case, the value of D_F(i) is set to the length of available following data minus CL or D_F(i)=t_ME−t_M-ACT(i)−CL, but not less than 0, where t_MEis the mapping-channel epoch length, in this example, 2000 msec, and CL is the reference-channel cycle length. Two such situations are illustrated in the example ofFIGS. 8A and 8B.
• The mapping-channel activation which is selected is the activation for which activation selection score A_SC(i) is a maximum. As shown inFIG. 8B,
• 
  A_SC(i)=t_NV(i)−D_P(i)−D_F(i).
• It is desirable that the selected mapping-channel activation be far in time from a ventricular-channel activation and that the neighboring cycle lengths in the mapping channel be close to reference-channel cycle length CL. This mathematical construction of the activation selection score A_SC(i) accomplishes this desired relationship.
• The computations outlined above and represented inFIG. 8B were carried out for the four candidate mapping-channel activations132bthrough132e(having activation times labeled130). The activations labeled132aand132fwere not detected;activation132aoccurred too close to the beginning ofepoch114, andactivation132fwas not large enough to trigger the threshold inactivation detection step74. Thus, as shownFIG. 8B, both the beginning-of-data rule and the end-of-data rule were applied in this example to determine values foractivations132band132e. In the row of data in the table for mapping-channel activation132b, the value of D_P=40 was determined by application of the beginning-of-data rule, and in the row of data in the table for mapping-channel activation132e, the value of D_F=247 was determined by application of the end-of-data rule.
• Mapping-channel activation132cis selected based on its maximum activation selection score A_SC=290 among the candidate mapping-channel activations.
• FIGS. 8C-1 through 8C-4 are a set of plots illustrating in detail method steps120,122, and124 in which LAT is computed based on selected mapping-channel activation132cand a reference-channel activation134. (Note that reference-channel activation134 is not detected in the method, but it is clear inFIG. 8A thatactivation134 is an activation near in time to t_M.) As indicated above, in this example, fiducial times t_Mand t_Rare the instants in the mapping-channel and reference-channel activations at which the maximum negative velocity occurs.FIG. 8C-1 illustrates an expanded signal of mapping-channel activation132c, andFIG. 8C-3 illustrates an expanded signal of mapping-channel activation velocity132c:v. Mapping-channel fiducial time t_Mis indicated by the dotted line labeled136. The value of t_Min this example is 4716 msec as indicated to the left of mapping-channel activation132cinFIG. 8A.
• Reference-channel activation134 is the activation in reference-channel signal108 which is located within ±CL/2 along the time axis of reference-channel signal108.FIG. 8C-2 illustrates an expanded signal of reference-channel activation134, andFIG. 8C-4 illustrates an expanded signal of reference-channel activation velocity134:v. Reference-channel fiducial time t_R, indicated by the dotted line labeled138, is clearly located within ±171 msec (±CL/2) of t_M.
• In this example, reference-channel activation134 occurs after mapping-channel activation132c, and the local activation time LAT=t_M−t_R=−15 msec. This value of LAT provides a single point in the generation of an LAT map. As mentioned above, an LAT map is based a single reference channel with its electrode placed at the same point in the cardiac structure throughout the entire generation of the map. A plurality of LAT measurements is used to generate an LAT map, each such point being made available for display by the system.
• InFIG. 7C, analternative embodiment122′ for LAT determination is illustrated.Alternative embodiment122′ takes advantage of the fact that multiple fiducial times t_Rare available in reference-channel signal108 to which mapping-channel fiducial time t_Mmay be compared. Inmethod step200, the times t_Rof maximum negative velocity generated inmethod step110 ofFIG. 7A are identified in reference-channel signal108. (As described above, activations in the reference channel are not detected using a threshold-crossing technique; rather a simple numerical search method may be used to find the times t_Rof local relative negative-velocity maxima insignal108.) Inmethod step202, the four nearest values t_Rto mapping-channel fiducial time t_Mare selected from among the values of t_R, and inmethod step204, each of these four times t_Rare adjusted by adding or subtracting multiples of reference-channel cycle length CL until each adjusted value t_RAis within ±CL/2 of t_Mso that t_RAnow represents its relative time position within one cycle length CL. Inmethod step206, these four values are averaged and this average av_RAis included with the set of four values of t_RAinmethod step208, creating a set of five time values. Inmethod step210, the median med_RAof this set is found and inmethod step212, LAT is determined as LAT=t_M−med_RA.
• Referring again toFIG. 8A, reference-channel signal108 is shown with a set ofnumerical values135 next to each activation insignal108. (Only onesuch triplet135 of values is labeled withreference number135 to avoid clutter withinFIG. 8A.) There are 14such triplets135 of values, and each triplet consists of (1) a fiducial time t_Rof the maximum negative velocity in signal108 (the number in parentheses), (2) the adjusted time t_RAby adding or subtracting multiples of cycle length CL to place t_RAwithin ±CL/2 of t_M, and (3) the time difference between t_Mand t_RA, also referred to as intermediate LAT values. This set of time values is also shown in the table ofFIG. 8D.
• Referring toFIG. 8D, rows A, B, and C correspond to the time value triplets described in the above paragraph. Row D is an ordered list of the intermediate LAT values in row C. Due to variations in measurement conditions or sources of variability, outlier values may be present at either end of this ordered list. To avoid such outliers (such as the two highest values in row D,44 and161), an interquartile set of values is selected in row E, dropping the lowest and highest 25% of values from the ordered list. This interquartile list may then be used to provide an estimate of a measurement confidence interval formethod step26 inFIG. 1. The range of this interval is indicated by the end values of row E which range from −15 to −11, or a ±2 msec LAT measurement confidence interval. Further, a measurement-confidence criterion formethod step26 inFIG. 1 may be as follows: if the width of the measurement confidence interval (in this example, 4) is greater than some percentage of cycle length CL, then alert the EP doctor inmethod step30. This width percentage criterion may be about 5% but this value is not intended to be limiting. An absolute width, say 15 msec, may also be used as the criteria; again such absolute width criterion value is not intended to be limiting.
• Referring toFIG. 8D to illustrate the alternative method of LAT determination described inFIG. 7C, the four times t_Rindicated byreference number214 are the four nearest times t_Rwhich encompass mapping-channel fiducial time t_M, forming aset214 of four values t_R(in row A) as described inmethod step202. The four values are adjusted as described above and form the set of four values (in row B) described inmethod step204. The average of these four values t is computed in method step206 (av_RA=4730.25 msec), This value av_RAis added to set214 to form a set of five values as indicated inmethod step208. The set of five values is now the set (4727, 4729, 4730.25, 4731, 4734). The median med_RAof this set is 4730.25 as found inmethod step210, and LAT is determined by LAT=t_M−med_RA, or LAT=4716−4730.25=−14.25 msec. One advantage of such an embodiment is that some additional computational precision occurs with the use of the averaging step.
• The use of the four nearest times t_Rwhich encompass t_Mis not intended to be limiting. Other choices for the number of values t_Rused in the LAT determination may be employed.
• Additionally, the steps described with respect toFIG. 8D also form a second alternative embodiment for LAT determination, potentially taking into account even more reference-channel fiducial times t_Rin the estimate of LAT. In this second alternative embodiment, the median value of the interquartile set of intermediate LAT values may be used as the LAT value for the current mapping point. As shown inFIG. 8D, the median value in row D is −13, the LAT value for this second alternative method of LAT determination. Again as above, other choices may be made for the number of values of t_Rused in the determination of LAT.
• Signal quality SQ as determined inmethod step79 ofFIG. 4A is also applicable for use withinmethod step24 ofFIG. 1 in which signal quality is monitored to provide alert28 when signal quality SQ degrades. One criterion by which to assess signal quality inmethod step79 ofFIG. 1 is simply to determine if all three signal quality values (ventricular channel, reference channel and mapping channel) are positive. The signal quality SQ determination for the mapping channel differs from that of the other two channels in that it does not include a signal irregularity SI term (i.e., SQ=SS−2N_Sfor the mapping channel) since the 2-second epoch of the mapping channel is too short to generate a meaningful value for signal irregularity SI. The decision ofmethod step79 is in the affirmative if any of the three signal quality SQ values is negative, at which time an alert is given to the EP doctor. Other criteria may be used inmethod step24 to triggeruser alert20.
• As described above, activation maps are used during certain cardiac procedures. But during such procedures, a variety of other cardiac parameters may advantageously displayed. Among these may be: (1) a value for starting reference-channel cycle length; (2) a value for current reference-channel cycle length CL with a confidence interval; (3) a value for LAT with a confidence interval; and (4) a value for ventricular-channel pulse interval PI with a confidence interval. All of these quantities are generated by the method disclosed herein. For example, a confidence interval for current reference-channel cycle length CL may be determined from the lag L2 of a peak in ACF near twice the cycle length CL, with the confidence interval being ±(L2−2CL) interval. A confidence interval for the LAT measurement may be ±half the interquartile range as described above. A confidence interval for ventricular pulse interval PI may be represented by range R_PI(±R_PI/2) as computed inmethod step92 ofFIG. 5.
• As described above, an activation map comprises a plurality of LAT measurements all of which are made relative to a particular reference-channel signal. One aspect of the inventive automatic method of measuring parameters of multi-channel cardiac electrogram signals includes the ability to compensate for signal degradation in the reference channel during the creation of an activation map. Since LAT maps are made relative to a specific reference channel, if the reference-channel signal being used degrades during mapping below a useful level of signal quality, the inventive method enables another reference channel to be selected and recreates the set of LAT measurements based on the new reference channel and generates a new map. This is possible since the inventive method computes reference-channel parameters as described above for several reference channels in real-time and stores the necessary parameters for use if needed. Very fast computation available with present computing equipment enables these “extra” channels to be recorded and analyzed in real-time without hindering the operation of the “current” channels being used to create a map.
• As seen above, a ventricular channel and a reference channel from among the channels of the MCCE signals are used in the automatic method of the present invention. The processes of selecting these channels automatically are among the various aspects of the inventive automatic method.FIG. 9 is a schematic diagram illustrating the inclusion of automatic selection of the ventricular and reference channels to the inventive automatic method of measuring parameters of MCCE signals. The steps of this overall combination are indicated by reference number140.
• Referring toFIG. 9, these automatic initialization steps include automatic selection of a ventricular channel, indicated byreference number144. The method steps of an embodiment ofautomatic process144 are illustrated inFIGS. 10A and 10B and described in detail below. Following the selection of a ventricular channel, a reference channel is automatically selected as indicated byreference number164, and the method steps of an embodiment ofautomatic process164 are illustrated inFIGS. 11A-11C.
• The entire automatic method of the invention disclosed herein is under the control of the electrophysiologist (EP doctor) as indicated above. At the time of a medical procedure, there may be overriding medical or technical reasons for the EP doctor to reject a channel or the channels which have been automatically selected, so automatic method140 includes aconfirmation step142 in which the EP doctor performing the procedure may accept or reject the channels which have been automatically selected. If the EP doctor rejects one or both of these selections, indicated by the “N” option onconfirmation step142, channel selection may be done manually or channels may be selected automatically as indicated bypathway142n.
• Upon final selection of ventricular and reference channels, automatic process140 continues with the method steps of mapping as indicated byreference number10 and as described in detail above.
• FIGS. 10A and 10B are schematic diagrams ofembodiment144 of the process of automatically selecting a ventricular channel from a set of candidate MCCE channels. The ventricular channel is selected from a set of candidate channels for use within the inventive automatic method of measuring parameters of MCCE signals. It is desirable that the ventricular channel selected be a channel which exhibits high signal quality and provides a stable representation of the action of the ventricles of the heart. As mentioned above, the ventricular channel of choice is most often connected to a body surface electrode, but other channels such as an epicardial electrode or an intracardiac electrode may provide the most useful signal among the set of candidate ventricular channels.
• Referring toFIG. 10A, in one embodiment ofautomatic process144 for selection of the ventricular channel, an initial time period for capturing and assessing signal characteristics from all candidate channels is employed. Data is captured during the time period represented by five epochs E1 through E5. (Epochs E1 through E5 are also called sub-signals E1 through E5.) In the embodiment ofFIG. 10A, this initial period is 30 seconds long. The signal from each candidate channel is an absolute-value velocity signal generated as described above. During the ventricular-channel selection process, there may be numerous possible channels being evaluated, the candidate set of channels including body-surface channels as well as some epicardial and intracardiac channels, depending on the strength of the ventricular signal in such channels. The process illustrated inFIG. 10A is applied to the waveform signals of each channel individually, generating for each channel two measures by which the ventricular-channel selection is made.
• Referring again toFIG. 10A, a 30-second signal148 from each candidate channel is divided into five six-second epochs (sub-signals) E1 through E5. The 30-second period, six-second epoch length and number of epochs in the initial period are not intended to be limiting.FIG. 10A illustrates the process for just one such candidate channel; a plurality of candidate channels is processed concurrently within ventricular-channel selection process144.
• As illustrated inFIG. 10A, the signals in each epoch E1-E5 are processed to determine the corresponding signal quality SQ_iinmethod step150. (Five such separate signal quality calculations are performed; only one is labeled withreference number150, but the marking SQ_iindicates this calculation being performed sequentially five times, each producing a value SQ_ifor the signal quality of the corresponding epoch E_i.)Signal quality calculation150 is carried out as illustrated inFIG. 4A.
• Inmethod step154, the five signal quality values SQ_iare summed to produce an overall signal quality value SQ_VCfor each candidate ventricular channel.
• Also illustrated inFIG. 10A, the signals in each epoch E1-E5 are processed in method steps152 to determine a value for pulse-interval variability VAR_ifor each of the five epochs (as found in the steps ofFIG. 5). Five such variability calculations are performed; one is labeled withreference number152, but the markings indicate that this calculation is performed sequentially five times, each producing a value VAR_ifor the variability of the corresponding epoch. In this embodiment, variability VAR_ifor each epoch is the difference between its maximum MAX_E1and minimum MIN_E1pulse-interval values during the six-second epoch. This is illustrated for epoch E1 in thedetailed method step152aas follows:
• 
  VAR₁=MAX_E1−MIN_E1
• Similar relationships for each epoch are calculated to generate the variability VAR_ifor each epoch E1 through E5.
• Inmethod step156, the maximum value of variability among the five values of variability is set as the variability VAR_VCof the candidate ventricular channel.
• At this stage in the automatic ventricular-channel selection process, each ventricular channel in the set of candidate ventricular channels has a channel signal quality assessment value SQ_VCand a channel pulse-interval variability assessment value VAR_VCwhich will be used to complete the automatic ventricular-channel selection process.
• FIG. 10B illustrates the final portion of ventricular-channel selection process144 for this embodiment. In general, it is desirable to select a ventricular channel with high signal quality and acceptable variability. Thus,FIG. 10B illustrates an embodiment which selects the channel with the highest signal quality but excludes channels from this choice which have especially high variability. Variability in such signals is expected. However, the reason that some channels are excluded is that these channels have variability values which are remarkably higher than those of other channels, indicating that the algorithm may be failing to detect activations properly in such channels.
• InFIG. 10B, the wide arrows labeled156a(two such arrows) and154erepresent a plurality of values as indicated.Arrows156atherefore represent that the channel variability values for each of the channels are all inputs to the method steps158 and160. Inmethod step158, a ventricular-channel variability threshold T_VCis calculated as follows:
• 
  T_VC=2·[median(VAR_VC)+Δ_VC]
• where Δ_VCis a small increment of time which may be added into this calculation simply as a computational convenience, such as to avoid singular calculations or to avoid excluding too many channels when the variability of some channels is extremely small. The inclusion of increment Δ_VCin the embodiment ofFIG. 10B is not intended be a limitation; Δ_VCmay have a value of 0 or other computationally-convenient value. InFIG. 10B, Δ_VChas a value of 5 milliseconds.
• Ventricular-channel variability threshold T_VCis a threshold value above which the variability of a channel is deemed to be unacceptably high. In method step160, the variability VAR_VCfor each channel is compared with ventricular-channel variability threshold T_VC, and channels for which VAR_VCis equal to or exceeds threshold T_VCare excluded from being the selected ventricular-channel VC_S.
• Other computational assessments of signal quality and variability for each channel and for the exclusion of channels on the basis of high variability are of course possible. The specifics of these assessment embodiments are not intended to be limiting.
• Wide arrow154erepresents one or more ventricular-channel signal quality values SQ_VCfor channels which have not been excluded in method step160. Each channel represented in the set ofvalues154eis a possible selected ventricular-channel VC_S. Inmethod step162, the channel with the highest value of channel signal quality SQ_VCis then selected as the ventricular-channel VC_Swithin the automatic method of measuring parameters of MCCE signals. (The method therefore also knows which channels are, for example, “second best” and “third best” among the candidate channels.)
• After ventricular-channel VC_Shas been selected using ECG signals over an initial period of time (30 seconds in the embodiment ofFIGS. 10A and 10B), MCCE signals from the initial period of time or from an additional period of time are used to select a reference channel to be used within the automatic method of measuring parameters of MCCE signals.FIGS. 11A through 11C are schematic diagrams ofembodiment164,166 of the process of automatically selecting a reference channel from a set of candidate channels.
• Body surface electrode channels are generally known not to be good choices for reference channels for many arrhythmias; thus, the reference channel is typically selected from the remaining set of MCCE channels for use within the automatic method of measuring parameters of MCCE signals. It is desirable that the reference channel selected be a channel which exhibits high signal quality and low cycle-length variability and also which exhibits a fast heart rate. For physiological reasons related to the cardiac measurements for which the present invention is intended to be used, it is also desirable that the selected reference channel indicate the shortest cycle length CL. All of these criteria are used to select a reference channel from among the set of candidate reference channels.
• Referring toFIG. 11A, inautomatic process164 for selection of the reference channel, an initial time period for capturing and assessing signal characteristics from all candidate channels is employed. Data are captured during the time period represented by five epochs, labeled E6 through E10. (Epochs E6 through E10 are also called sub-signals E6 through E10.) In the embodiment ofFIG. 11A, this additional initial period is 30 seconds long. The signal from each candidate reference channel is an absolute-value velocity signal generated as described above. During the reference-channel selection process, there may be numerous possible channels being evaluated. What specific channels are candidate reference channels for cardiac mapping are well-known to those skilled in electrophysiology. The process illustrated inFIG. 11A is applied to the waveform signals of each candidate channel individually, generating for each channel two measures by which the reference-channel selection is made.
• Referring again toFIG. 11A, a 30-second signal168 from each candidate channel is divided into five six-second epochs (sub-signals) E6 through E10. The 30-second period, six-second epoch length and number of epochs in the initial period are not intended to be limiting in any way.FIG. 11A illustrates the process for just one such candidate channel; a plurality of candidate channels is processed concurrently within reference-channel selection process164.
• As illustrated inFIG. 11A, the signals in each epoch E6-E10 are processed to determine the corresponding signal quality SQ_iinmethod step170. (Five such separate signal quality calculations are performed; only one is labeled withreference number170, but the marking SQ_iindicates this calculation being performed sequentially five times, each producing a value SQ_ifor the signal quality of the corresponding epoch E_i.)Signal quality calculation170 is carried out as illustrated in the steps ofFIG. 4A.
• Inmethod step174, the five signal quality values SQ_iare summed to produce an overall signal quality value SQ_RCfor each candidate reference channel.
• Also illustrated inFIG. 11A, the signals in each epoch E6-E10 are processed in method steps172 to determine a cycle length CL_iand a value VAR_ifor cycle-length variability for each of the five epochs. The cycle length determination is performed using the method steps described in detail above and illustrated inFIGS. 6A and 6B (and as referred to inFIG. 11B). The determination of cycle-length variability (and cycle length) for each epoch is described using the schematic diagram ofFIG. 11B, which illustrates thisprocess172afor epoch E6 as indicated by the dotted-line ellipse inFIG. 11A. Five such cycle length and cycle-length variability calculations as illustrated inFIG. 11B are performed; one is labeled withreference number172a, but the markings of each such element inFIG. 11A indicate that this calculation is performed sequentially five times, each producing a value VAR_ifor the cycle-length variability of the corresponding epoch.
• In the method steps illustrated inFIG. 11B, F6 is the magnitude-coincidence autocorrelation function of epoch signal E6 as calculated inmethod step178. The threshold value for the magnitude-coincidence autocorrelation is a function of the noise in signal being processed as described in the summary section which defines magnitude-coincidence autocorrelation and with respect tomethod step226 ofFIG. 6A. Cycle length CL is computed inmethod step182 as indicated, using steps outlined inFIGS. 6A and 6B. Then, inmethod step186, a peak in F6 near a value of lag near 2 times CL is identified, and an interim variable DCL is set to the lag at such peak. Inmethod step188, another interim variable CL_Ais set to DCL−CL. The value for the variability VAR₆of epoch E6 is then calculated as follows in method step190:
• 
  VAR₆=|CL−CL_A|
• Similar relationships for each epoch are calculated to generate the variability VAR_ifor each epoch E6 through E10.
• Referring again toFIG. 11A, inmethod step176, three additional values (in addition to SQ_RC) are determined for each epoch E6 through E10 by which to select from among the candidate reference channels: VAR_RC; MAX_CL; and MIN_CL. MAX_CLis the maximum cycle length CL_iamong the five cycle-length values. MIN_CLis the minimum cycle length CL_iamong the five cycle-length values. VAR_RCis the maximum variability value VAR_iamong the five variability values.
• At this stage in the automatic reference-channel selection process, each reference channel in the set of candidate reference channels has a channel signal quality assessment value SQ_RC, a channel variability assessment value VAR_RC, and maximum and minimum cycle length values MAX_CLand MIN_CLwhich will be used to complete the automatic reference-channel selection process.
• FIG. 11C illustrates the final portion of reference-channel selection process164 for this embodiment. In general, it is desirable to select a reference channel with high signal quality, low variability, and short cycle length, as indicated above.FIG. 11C illustrates an embodiment which selects the channel with these characteristics.
• InFIG. 11C, the wide arrows labeled174a,176v,176x, and176neach represent a plurality of values as indicated.Wide arrow174arepresents all values of SQ_RCfrom the candidate reference channels;wide arrow176vrepresents all values of VAR_RCfrom the candidate reference channels;wide arrow176xrepresents all values of MAX_CLfrom the candidate reference channels; andwide arrow176nrepresents all values of MIN_CLfrom the candidate reference channels.
• Inmethod step192, a figure-of-merit FM_RCis evaluated for each candidate reference channel. FM_RCfor each candidate reference channel is computed as follows:
• 
  FM_RC=SQ_RC/S_RC−MAX_RC−MIN_RC−S_VAR·VAR_RC
• where S_RCis an arbitrary scale factor and S_VARis an arbitrary scale factor. The two scale factors are chosen such that a useful tradeoff within the figure-of-merit FM_RCis created. When signal quality values SQ_RCare in microvolts and cycle lengths are in milliseconds, a value of S_RCof 32 and a value of S_VARof 2 have been found to yield a useful tradeoff among cycle lengths, variability, and signal quality and also to be computationally convenient.
• The FM_RCvalues for each candidate reference channel are output frommethod step192 as indicated bywide arrow194. Inmethod step196, the channel with the highest value of FM_RCis the selected reference-channel RC_S.
• Other computational assessments of signal quality, variability, and cycle length for each channel are of course possible. The specifics of these assessment embodiments are not intended to be limiting.
• As described above, one aspect of the inventive automatic method of measuring parameters of multi-channel cardiac electrogram signals includes the ability to compensate for signal degradation in the reference channel during the creation of an activation map by selecting a new reference channel and recreating the set of LAT measurements based on the new reference channel and generating a new map. During the initial selection process for the reference channel, the inventive method keeps track of the reference channels which have values for FM_RCjust below the selected reference channel RC_Sso that if necessary, these “second best” reference channels can be substituted for the selected reference channel and the mapping process can continue without losing the valuable time and effort that has already been spent on the mapping process.
• In another aspect of the inventive method, multiple mapping channels may also be employed, and the processing steps outlined herein applied to multiple mapping channels as well as multiple reference channels. Some catheters which are used in cardiac procedures may include multiple electrodes in a variety of configurations. In addition, multiple catheters may be employed. The speed of computer processing available enables numerous calculations to be done very rapidly such that multiple mapping channels may be supported to generate a plurality of maps as the EP doctor moves the mapping electrodes throughout chambers and vessels of the heart.
• FIG. 12 is a matrix which schematically illustrates a series of reference channels and mapping channels from among a set of MCCE signals which in an aspect of the inventive method may be processed in parallel to generate multiple LAT maps as well as track the other cardiac parameters measured with the inventive method in a variety of combinations of reference and mapping channels. Some channels may be displayed for the EP doctor in parallel while others may serve as possible backups in case of signal degradation as explained above.
• Referring toFIG. 12, the table260 shown is an 8×8 matrix (eight reference channels R₁through R₈and eight mapping channels M₁through M₈. The number of channels illustrated is not intended to be limiting, but rather to illustrate the concept of flexible configuration combinations of reference and mapping channels. In fact, a set of MCCE signals is typically much larger than the eight illustrated here. Eight channels in the MCCE set of signals is only for purposes of illustration.) The “

  

  ” symbols indicate that the mapping channel and reference channel at such an intersection together can generate a map. In other words, as shown in table260 ofFIG. 12, each channel M₁through M₈can be paired with any other channel R₁through R₈, and any combination of these pairings can be created. The “-” simply indicate redundancy in table260. The “-” symbol is located along the diagonal of table260 indicating that any channel is not usefully paired with itself. Further, some channels among the set of channels may in fact be ventricular channels (primarily exhibiting ventricular signal characteristics) and therefore also not good candidates for either a reference or mapping channel.
• The advantages of such multiple-channel processing configurations are that procedure time may be shortened but also that a much richer array of measurements may be obtained to provide better information to the EP doctor to ameliorate the cardiac deficiency being treated. Further, as described above, backup channels can be available to deal with lost or degraded signals during a procedure without the need to start the procedure over again.
• It is possible in some multi-channel configurations that certain information may be shared among several parallel computations. For example, it is possible that ventricular pulse-interval values may be used for the determination of several reference-channel cycle lengths, and ventricular-channel activation times may be shared for use with more than one mapping channel. And many other combinations other than those exampled here are possible with the multi-channel processing of the inventive method described herein.
• Inembodiment10 of the inventive automatic method of measuring parameters of multichannel cardiac signals described in detail above, contiguous six-second epochs of MCCE signal data are used. Alternatively, a moving-window format of selecting epoch starting and end points may be used, such as the next epoch in a series of epochs consisting of the last 5 seconds of the previous epoch and a new sixth second. Other epoch-forming strategies may be used, depending on the computational advantages which may be possible and on the desired speed of analysis.
• FIG. 13 is a schematic block diagram illustrating an alternative embodiment of the monitoring of cardiac channel quality. The alternative embodiment replaces a portion of the schematic block diagram ofFIG. 1. As previously described and as illustrated inFIG. 1 and inFIG. 4A, signal quality is defined as a quantity dependent on signal strength, signal irregularity and noise. In addition to such a specific assessment of signal quality, an overall channel-quality assessment may include information on the variability of the cycle length (or pulse interval when considering ventricular channels). Such an overall assessment for a reference is represented by a figure-of-merit which reflects the fact that higher signal quality values and lower cycle-length variability are desirable.
• In order to make an assessment of overall channel quality during ongoing operation of the inventive method, one embodiment of such a system includes applying the steps of automatic channel selection for initialization in real-time to monitor channel quality. These channel selection steps are fully described above and illustrated inFIGS. 9-11C. As described in these embodiments, assessment of overall channel quality includes calculation of various signal quantities over a plurality of epochs. In these embodiments, five contiguous six-second epochs are evaluated as shown inFIGS. 10A-11C. The number of epochs may vary and the length of the epochs may also differ from 6 seconds. Further, signal parameters may be computed over moving time windows rather than sequential (contiguous) epochs. The epochs may be of fixed-length or of various lengths. Any number of variations of time periods are possible in the computation of the elements of channel figures-of-merit, signal quality and cycle-length or pulse-interval variability.
• The schematic ofFIG. 13, labeled withreference number24a, replaceselement24 inFIG. 1 in the alternative embodiment ofFIG. 13. The connections to the block diagram are shown as connections withelements22,26 and28 ofFIG. 1. Thefunctional elements25 and27 represent the alternative functions of calculating a figure-of-merit FM (element25) and using figure-of-merit FM in comparison to a channel quality standard (element27). Such a standard may be a preset threshold value that figure-of-merit FM must exceed or it may be a dynamically-determined standard such as being the highest figure-of-merit FM among a set of candidate channels. Other standards are also possible such as, for example, a combination of a preset threshold value and the maximum among a set of channels. In such a case, a channel has acceptable overall quality when it has a figure-of-merit FM which is both the highest in a set of channels and is above a preset threshold value.
• The calculation of figure-of-merit FM as illustrated inFIG. 13 involves the application of the steps ofFIGS. 10A and 10B for the assessment of overall quality for candidate ventricular channels and application of the steps ofFIGS. 11A-11C for the assessment of the overall quality for candidate reference channels. As described above, various strategies for determining the timing of signals used are possible.
• In this inventive method, a plurality of channels are stored and processed such that the monitoring of overall channel quality is possible for cardiac channels as desired. The inventive method includes a variety of strategies for monitoring overall cardiac channel quality, including (a) performing the necessary calculations in real-time for only one or more of the “active” channels (the currently used mapping, ventricular, and reference channels for LAT determination), (b) performing the necessary calculations in real-time for the entire plurality of cardiac channels, and (c) performing such calculations for a subset of cardiac channels. Since channel signal data is stored, if a strategy such as (a) is chosen, overall channel quality of other (“non-active”) channels can be performed when necessary to determine which channel will replace the current cardiac channel. If a strategy such as (c) is employed, the inventive method monitors every cardiac channel, making an up-to-date assessment of the overall channel quality for every cardiac channel in the system available at any time.
• FIG. 14 is a schematic block diagram illustrating an alternative embodiment of the channel selection method illustrated inFIG. 9, adding to real-time operation the steps of automatic channel selection which in the embodiment ofFIG. 9 operate only within the initialization process. The method steps ofFIG. 14, indicated together as140a, include not only the initialization steps from the embodiment ofFIG. 9 but alsoelements202 through210.Element10a, as shown, is a modification ofelement10 inFIG. 9, modified to include the now real-time automatic steps embodied inFIGS. 10A-11C and 13.Elements202 through210 describe an embodiment which carries out cardiac channel replacement during real-time operation of the inventive method, replacing a reference or ventricular channel when the overall quality of such channel degrades during operation.
• Method element10aencompasses the running (real-time) operation of the inventive method as detailed inFIGS. 1-11C and 13.Connection202 andmethod element204 illustrate the step of deciding if a cardiac channel needs to be replaced, the result of calculating figure-of-merit FM instep25 and making a comparison with a channel-quality standard in step27 (seeFIG. 13). A user alert is generated in step28 (seeFIG. 1), anddecision step204 may proceed in a completely automatic fashion if the user chooses such automatic functioning. In such a case, an alert atstep28 is still provided to the user. However, it is anticipated that the user may intervene in the automatic steps of channel replacement to exercise more control over the replacement process. An “N” decision indecision step204 returns the process to overall running-time operation represented inelement10a.
• The inventive method inembodiment140aproceeds tomethod element206 in which a replacement channel is selected, either automatically or manually by the user, based on assessments of possible replacement cardiac channels.Embodiment140aproceeds then tomethod step208 in which confirmation of the channel replacement is carried out, either automatically or by user intervention. With a “Y” decision atdecision step208, the process continues channel replacement and updating inelement210 and then the process continues with overall running-time operation represented inelement10a. An “N” decision atdecision step208 returns the process to overall running-time operation represented inelement10awithout channel replacement and updating.
• In the embodiments ofFIGS. 10A through 11C, the assessments of overall channel signal quality differ for candidate reference channels and candidate ventricular channels. Skilled practitioners in the use of cardiac channel signals and in electrophysiology are able to determine which cardiac channels may serve as possible reference channels as well as those cardiac channels which have the properties of ventricular channels. Herein, the set of candidate reference channels is also referred to as a first subset of cardiac channels, and the set of candidate ventricular channels is also referred to as a second subset of cardiac channels.
Transformation of Local Activation Times
• FIG. 15 is a high-level schematic block diagram illustrating the steps of anembodiment300 of the inventive method for transforming LAT values using a second reference channel signal when the timing stability of a first reference channel signal degrades below a timing-stability standard, in order to avoid substantial loss of LAT values in spite of the loss of timing stability.Method element302 represents a plurality of reference channel signals among the MCCE electrogram signals available for processing. The ventricular channel and mapping channel signals306 are also among this set of available MCCE signals.
• Method elements304 (three shown) represent the step of storing these cardiac signals.Method embodiment300 may be realized within a computer programmed to carry out the steps as described herein, and the divisions between the method elements may vary based on the programming generated to carry out such steps. Other quantities in addition to the cardiac electrogram signals may be captured in computer memory, although not shown inFIG. 15 (other than timing offsets in method element308), in order to enable other computed values to be generated as desired, either during a procedure or retrospectively.
• Connection302cis drawn as a wider arrow to represent that the plurality ofreference channels302 are available within the method structure along aconnection302cas shown. A connection labeled302cis shown in three other locations within the schematic block diagram ofFIG. 15 to indicate thatreference channels302 are also available at these points within the block diagram ofembodiment300. In similar fashion, other connections are indicated by reference numbers ending with the letter “c” to indicate that multiple quantities are available thereon as indicated by the reference number of the method element and/or the nature of the multiple quantities available at such connections. For example,method element314, to be described later, is the source for the quantities available onconnection314c.
• Three method elements inFIGS. 15 (310,314 and318) involve the application of a standard by which a determination is made withinembodiment300 of the inventive automatic LAT-determining method. These standards are indicated by reference numbers ending with the letter “s” and which include the reference number of the method element in which such a determination is made.
• A first reference channel is established as the base reference channel for the LAT computation carried out inmethod element320.Method element320 is shown as having the three requisite channel signals (ventricular, mapping and reference) as indicated by connections frommethod elements306 and318. LAT computations are carried out inmethod element320 using a number of possible ways of determining parameters within such three signals. One such approach to LAT computation is described in this document inFIGS. 7A and 7B, their related figures, and the accompanying description. However, this example approach is not intended to be limiting; numerous other approaches for determining the necessary quantities for an LAT computation are within the scope of the present invention.
• Inmethod element308, timing offsets are computed to be used in the determination of LAT values as required and also for a determination inmethod decision element310 of whether or not the timing stability of the first reference channel signal (base reference channel signal) has fallen below a timing-stability standard310s. (Timing-stability standard310swill be described later in this document.) Inembodiment300, the timing offsets used to transform local activation times are themselves also local activation times and are defined as LAT_K(J) where LAT_K(J) is the local activation time of a reference channel J based on a reference channel K. Inmethod element308, timing offsets are computed for pairs of the plurality ofreference channels302. Aconnection306vis so labeled to indicate that only the ventricular channel signal is needed for the computations of timing offsets inmethod element308.
• For example, if there are fivereference channels302, then there are ten pairs withinreference channels302 for which timing offsets may be calculated. (In general, with n possible reference channels available, there are (n²−n)/2 timing offsets which may be computed.) Note that LAT_J(K)=−LAT_K(J) so that in this case there are 20 timing offsets which may be available for later use if each pair has had a period of time during which the timing stability of both signals has been satisfactory. It is most desirable that all possible pairs may be used for such timing offsets computations, but alternatively only a subset of the available pairs may be used.
• Inmethod decision element310, if it is determined that timing-stability standard310shas been met,method embodiment300 simply proceeds along the “Y” branch ofmethod decision element310 to process the next epoch of electrogram signal data inmethod element320 without replacing the first reference channel (the base reference channel) with a second reference channel. However, if it is determined that timing-stability standard310shas not been met (i.e., there likely has been a disruption caused, for example, by electrode displacement, degraded electrode contact or noise), the method ofembodiment300 proceeds along the “N” branch ofmethod decision element310 tomethod element312 in which certain characteristics of the reference channel signals302 are computed to enableembodiment300 to make a determination of which channels among the reference channel signals302 have been disrupted.
• Method decision element310 as described here, while determining that timing stability has been degraded below timing-stability standard310s, does not determine which channel(s) has/have caused this loss of timing stability. Such determination is made inmethod elements312 and314. Note that four specific embodiments of signalcharacteristic computations312 are described later in this document. These four embodiments are not intended to be limiting; other computations of signal characteristics may be used within the intent ofmethod element312.Connection312cindicates that the signal characteristic results ofmethod element312 are available tomethod element314.
• Method element314 represent the method step of sortingreference channels302 according to whether or not specific channels have been disrupted and are therefore not available to be selected inmethod step318 as the second reference channel. This sorting is accomplished by comparison of computed signal characteristics with a channel-sorting standard314s. The output ofmethod element314 is subset of reference channels available alongconnection314cwhich may be selected as the second reference channel (new base reference channel for LAT determination).
• Inmethod step316, the signal qualities of each available reference channel are determined, and inmethod step318, a selection is made of the second reference channel based on a signal-quality standard318s. A number of signal quality computations may be made inmethod element316, examples of which are high signal amplitude, low signal noise, low signal amplitude variability, low cycle-length variability and short cycle length. Some approaches to computing signal quality are described earlier in this document, such as in sections related to automatic selection of channel signals. None of these specific measures of signal quality are intended to be limiting to the scope of this invention; other computations of signal quality are possible within the intent ofmethod element316.Connection316cindicates that the signal characteristic results ofmethod element312 are available tomethod element318.
• During a cardiac mapping procedure, the ventricular channel electrode and the electrodes of thereference channels302 ideally remain stationary. Since the reference-channel electrodes ideally remain stationary, a set of timing offsets LAT_J(K) are computed inmethod element308 using pairs from among thereference channels302, the values of the timing offsets being computed in advance of their being needed.
• MCCE signals are by their very nature noisy from a number of sources, and thus the determination of LAT values, and of course, timing offsets, is a statistical process. Transformation factors may be established by averaging a set of LAT_J(K) values for the pair of reference channels (J,K) over a period of time during a cardiac mapping procedure, and the statistical variations in LAT_J(K) may be used as part of a determination as to whether a reference channel is suitable to remain as the first (base) reference channel or to be available to be selected as the second (new base) reference channel. Among the causes of stability degradation is of course the physical movement of a reference channel electrode such that its location has changed and therefore cannot be relied upon as a suitable reference.
• When a new base channel (second reference channel) has been selected inmethod element318, the computation of LAT values inmethod element320 changes from the straightforward use of the ventricular channel, the mapping channel and the base reference channel to include the use of one or more timing offsets. The availability of the timing offsets inmethod element320 is indicated by the timing offsets being inputs alongconnection308c.
• The transformation process triggered by a change from the first reference channel to the second reference channel is represented by the equation LAT_A(C)=LAT_A(B)+LAT_B(C) where timing offset LAT_B(C) is the local activation time of reference channel C based on a reference channel B, and such timing offset LAT_B(C) transforms values LAT_A(B) into LAT values LAT_A(C). As indicated, computation of LAT values inmethod element320 may include a straightforward calculation of an LAT value without the application of a timing offset or it may include transforming a past or current LAT value with the appropriate timing offset.Method element322 represents that LAT values may be displayed in a map or other form as desired by the EP doctor.
• One important aspect of the inventive automatic method is that it is able to transform both past and future LAT values so that an LAT map, once it has been started and considerable investment of time (both patient and medical personnel) has been undertaken along with possible exposure to X-ray radiation, the LAT map is able to be completed without the loss of such investment made in the capturing of LAT values.
• As an example, an LAT map is partially generated based on the first reference signal, such LAT values of mapping channel M being represented by LAT₁(M). Then, during this procedure the quality of the first reference signal degrades, and the first reference channel signal is replaced by a second reference channel signal. Each of the past LAT values LAT₁(M) may be transformed into LAT values LAT₂(M) by LAT₂(M)=LAT₂(1)+LAT₁(M). New LAT values LAT₂(M) may continue to be captured based on the second reference channel signal, and the completed map in this case is now based on the second reference channel. However, all of the points LAT₁(M) that were captured based on the first reference channel have been used to build the new map without having to be recaptured. The timing offset LAT₂(1) has been established over time during periods in which both the first and second reference channel signals were signals having suitable timing stability and quality.
• In another example, an LAT map is partially generated based on the first reference signal, such LAT values of mapping channel M being represented by LAT₁(M). Then, during this procedure the quality of the first reference signal degrades, and the first reference channel signal is replaced by a second reference channel signal. New LAT values LAT₁(M) may be determined by computing LAT values LAT₂(M) based on a second reference signal and then transformed into LAT values LAT₁(M) by LAT₁(M)=LAT₁(2)+LAT₂(M). The completed map in this case is still based on the first reference channel, and all of the points LAT₁(M) that were captured based on the first reference channel have been used in the map with new LAT values added to the map having been transformed into LAT values based on the same (first) reference signal. The timing offset LAT₁(2) has been established over time during periods in which both the first and second reference channel signals were signals having suitable timing stability and quality.
• In the high-level schematic block diagram ofFIG. 15, a number of more detailed steps are implicitly embedded in the method elements shown inembodiment300. For example, among these functions are: (1) the logic and control withinmethod element320 to determine whether a future-value transformation or a past-value transformation is appropriate for the LAT computation; (2) the sequential processing of the channel signals (e.g., from epoch to epoch); and (3) the logic and control necessary to capture an LAT value when the EP doctor decides that such a value should be captured. The automatic nature of the inventive method of determining LAT is not compromised by the manual intervention of the EP doctor to “take a point” for mapping. The method proceeds automatically, computing LAT values, monitoring timing stability, and selecting replacement channels as described herein, with a “take point” process operating over the top of such automatic method.
• Also implicitly included in the method elements ofembodiment300 is logic control of the automatic method whereby after a reference channel has been replaced due to a loss of timing stability, a period of time is provided during which the statistics of the timing offsets are reestablished with the new set of electrode conditions such as might have been caused by a reference-channel electrode having shifted its position. For example, if six-second epochs are being used and twenty data points is regarded as sufficient, a period of two minutes would be required to reestablish such statistics.
• FIG. 16 is a table309 showing exemplary mean values for timing offsets as computed inmethod element308 ofFIG. 15. This example is of a reduced set of six cardiac channels. As described above, the timing offsets are local activation times LAT_J(K), the local activation time of a reference channel K based on a reference channel J where channels J and K are reference channels within a set of possible or candidate reference channels in the reduced set of MCCE channels. In this example, none of these channels would intentionally be used as an exploring catheter (mapping channel).
• The first two columns from the left of table309 indicate the channel numbers for the timing offset calculations. Each entry in the third column of table309 contains an LAT value to be used as a timing offset for the activation delay from one possible reference electrode to another possible reference electrode. As such, reference channel K takes the role of the mapping channel within an LAT calculation inmethod element308, and a negative LAT value indicates that the activation arrives first at the channel K electrode which is filling the role of mapping channel. As mentioned above, LAT_J(K)=−LAT_K(J); thus, half of the available timing offsets need not be shown in table309. Since six possible reference channels are available in this example, fifteen timing offsets are computed, and thirty offsets are available.
• The timing offset values shown in the third column in table309 are the mean values of each computed timing offset, averaged over a predetermined number of computed values. For this example, the timing offsets are computed for six-second epochs of signal data over a period of two minutes; therefore, 20 timing offset values are averaged to produce the means values shown in table309. The predetermined number (20) of epochs and the 6-second length for each epoch are not intended to be limiting; other values for these two variables may be used.
• It is assumed that repeated measurements of the timing offset between two channels with stationary electrodes will populate a Gaussian distribution with a mean and standard deviation (stdev). The fourth column of table309 shows the standard deviation of each of the 15 timing offset values over the 20-sample (two-minute) history of computed timing offset values. The statistics of the third and fourth columns are determined on a two-minute moving window to be used in the determination ofmethod element310 as to whether or not timing stability has been lost. Data from the epoch immediately following the two-minute window are used in conjunction with the statistics of the two-minute moving window to makedetermination310. If it is determined inmethod element310 that timing stability has not been lost, the means and standard deviations for the new two-minute moving window are updated and are ready for the data from the next six-second epoch to be processed.
• The data in the fifth and sixth columns represent a next epoch (epoch following the updating of the moving-window statistics). These two columns represent data taken from an epoch in which timing stability has not been lost, i.e., before loss of timing stability. The fifth column of data, labeled X, contains 15 individual computed timing offset values, one for each pair of channels in the set of timing offset pairs. For example, the timing offset LAT₁(2)=−24 msec and the timing offset LAT₃(6)=−42 msec.
• Method element310sindicates the use of a timing-stability standard (also indicated byreference number310s) in the determination ofmethod element310. Since it is assumed that the timing offsets have Gaussian distributions, one example of a useful timing-stability standard310sis to compute Z-scores for the timing offsets X for the epoch immediately after the current two-minute moving window for which the statistics have been calculated. This statistical test quantity is computed as
• 
  Z-score=(X−mean)/stdev
• and values for the fifth column of timing offsets X are shown in the sixth column. A Z-score exceeding ±2 indicates that a new timing offset (in the set of values labeled X) has significantly (with 95% confidence) deviated from the expected distribution, thereby indicating that some change has occurred. As can be seen in the sixth column, none of the Z-scores exceeds ±2, indicating that it is unlikely that timing stability has been lost as defined by this example timing-stability standard310s.
• The seventh and eighth columns of table309 are outlined in bold. The data in the seventh column of table309 illustrates an example of an epoch before which timing stability has likely been lost. Asmethod embodiment300 proceeds and an epoch has been found not to have indicated a loss of timing stability, the two-minute moving-window statistics would be updated with the latest data before proceeding to the next epoch data. In the simplified example of table309, the timing offsets values X in the seventh columns are assessed against the statistics in columns three and four to illustrate the detection of the loss of timing stability.
• The eighth column shows the Z-scores for the timing offset values in the seventh column, and it can be seen that two timing offset values, LAT₂(4) and LAT₂(6), have Z-scores which exceed −2, indicating that timing stability is very likely to have been lost prior to this epoch. (These two Z-scores are highlighted by gray shading.) Because timing offsets represent a relationship between two channel signals, when a timing offset changes, it is not known which channel (or both) have shifted or otherwise have lost a satisfactory signal. Inmethod element312 ofembodiment300, signal characteristics of reference channel signals are computed for the purpose of determining which channel(s) has/have caused the loss of timing stability as detected inmethod element310.
• When timing stability has been found to have been lost, it may become necessary to provide extra time during which the statistics of some computed timing offsets are reestablished. This could be achieved by enabling a timer to prevent such channels from use during, for example, a two-minute period.
• The timing stability standard310sexampled in table309 is not meant to be limiting; other indications of loss of timing stability may be used and are within the scope of the claimed invention.
• FIG. 17 is schematic block diagram illustrating the steps of an FFT-based method embodiment312(1) for computation of a signal characteristic for use withinmethod embodiment300. This signal characteristic computation312(1) generates an FFT-based parameter called the epoch center-of-power frequency freq_COPand is used within the example ofFIGS. 16 through17N illustrating embodiment300. In the method steps of embodiment312(1), a fast Fourier transform (FFT) is computed for each relevant channel signal in an epoch. Computation of an FFT is well-known to those skilled in the art of signal processing and is not detailed herein.
• FIG. 17A is an exemplary six-second epoch of a representative cardiac channel electrogram signal which is used to illustratemethod embodiment300, and in particular, to illustrate four alternative methods312(1)-312(4) for computation of a signal characteristic for use withinmethod embodiment300. (The abscissa for the exemplary epoch ofFIG. 17A is in millisecond (msec) or sample index for the 6,000 signal values in a six-second epoch.) Inmethod element324 ofFIG. 17, each such channel signal is divided into overlapping segments of data for the purpose of a more time-efficient computation of the FFT—by computing an FFT for each of a plurality of shorter overlapping segments which span the time range of the signal. (This method for estimating the power of a signal as a function of frequency is similar to Welch's method or Welch's periodogram.)
• A six-second window of signal data with a sampling rate of 1,000 samples per second (sps) is divided into five overlapping blocks of 2,048 values of the signal which span the entire 6,000 data values in the six-second epoch, as illustrated inFIG. 17B. Each of the five curves326(1)-326(5) is a raised and scaled sinusoidal curve also illustrating the weightings applied inmethod step326. The weightings326(1)-326(5) are applied individually to segments of the signal data. Each curve has 2,048 non-zero values. InFIG. 17B, curve326(1) spans betweendata values 1 and 2,048; curve326(2) spans between data values 977 and 3,024; curve326(3) spans between data values 1,977 and 4,024; curve326(4) spans between data values 2,977 and 5,024; and curve326(5) spans between data values 3,953 and 6,000.
• FIGS. 17C-17G are plots of the resulting weighted segments of signal data. (Note that the abscissa of each of the five plots inFIGS. 17C-17G span from 1 through 2,048 since the fast Fourier transform to be performed is on each segment individually although the weighted signal data is taken from different segments of the epoch as indicated above.) The segmenting and the application of weightings as illustrated in this example causes a net attenuation of the total power in each segment since the maximum weighting in each segment is 1.0 and the average across each segment is less than 1.0. However, since the same computation is applied to each epoch, the relative relationship of the resulting values is preserved, and the computation is much more time efficient.
• Inmethod element328, an FFT is computed for each of the weighted signal segments. Each of the five segment FFT's computed inmethod element328 is a 1,024-point FFT, i.e., there are 1,024 frequency coefficients (bins) each being 0.488 Hz wide (1,000/2,048). With a sampling rate of 1,000 sps, frequencies of up to 500 Hz can be measured. The resulting FFT coefficients are complex numbers, and inmethod element330, these complex coefficients are converted to magnitudes to form 1,024 real coefficients.
• FIGS. 17H-17L show the five resulting segment spectra for the weighted segment signal data ofFIGS. 17C-17G, respectively, including the conversion of the FFT coefficients to magnitudes inmethod element330. The average of these individual spectra is computed inmethod element332 and results in the spectrum ofFIG. 17M for the entire six-second epoch of signal data. (For computational efficiency, a sum could also be computed inmethod element332 instead of an average, the only difference being the scale factor; these approaches are equivalent.)
• Inmethod element334, the first moment of power is computed in order to determine the center-of-power frequency (freq_COP), as follows:
• 
  freq_COP=0.488·[Σ(c_i²·i]/Σc_i²
• where the c_i's equal the 1,025 spectral coefficients (the 0 Hz and 0.488 Hz coefficients are set to 0 to remove these frequencies), i equals the bin number of the frequency bin in the spectrum, ranging from 0 to 1,024, and the sums are computed over the range of 1,025 bins (the abscissa of each plot is bin number). (i=1 is the 0.488 Hz bin number.) Conceptually, this center-of-power frequency can be thought of as the average of the bin distances (bin number) from 0 Hz with the bin distances being weighted by the power at the corresponding frequency of the bin. The corresponding frequency of each bin is 0.488 times the bin number so the maximum frequency available in the 1,000 sps signal is 500 Hz (the Nyquist frequency). Note that the 0 Hz and 0.488 Hz (bin numbers 0 and 1) coefficients are set to zero since these frequencies are substantially adulterated by the use of the raised-cosine segmentation process and therefore are not useful for determining the power spectrum.
• The value of freq_COPfound for the spectrum ofFIG. 17M is 107 Hz. Note that the exemplary signal epoch data ofFIG. 17A has been captured through a 60 Hz notch filter as is often the case for cardiac electrogram signals. This can be seen inFIGS. 17H-17M as indicated by the deep V cut into the spectra around 60 Hz. As long as the filter remains in place throughout a procedure, its use is of little consequence to the ability of the inventive method to discriminate between signals.
• It has been determined that center-of-power-frequency freq_COPprovides an excellent indication that a change has occurred in a cardiac channel signal.FIG. 17N shows a table311 of the results of computations inmethod element312 of the FFT-based freq_COPsignal characteristic for the example reduced channel set ofFIG. 16. In the first column of table311, further specifics of the six channels are given.Cardiac electrogram channel 1 is shown as HRA (high right atrium) andchannel 2 as HIS (bundle of His). These two channels are independent of all of the cardiac channels. (The bundle of His is a slender bundle of modified cardiac muscle that passes from the atrioventricular node in the right atrium to the right and left ventricles by way of the septum and that maintains the normal sequence of the heartbeat by conducting the wave of excitation from the right atrium to the ventricles. It is also called the atrioventricular bundle.)
• The remaining four channels (A12, A34, A56, and A78) are all bipolar signals between two adjacent electrodes of an octo-polar catheter. (The signal epoch data ofFIG. 17A is data from a representative six-second epoch from channel 5-A56.) The electrodes are at fixed spacings on the catheter and therefore are generally constrained to move together. The second column of table311, labeled “Wire,” shows an index of wires. Channels having the same wire index are on the same wire or catheter, and all channels having the same wire index are likely to be disrupted if the wire or catheter is disrupted.
• The third and fourth columns of table311 present the statistics (mean and standard deviation stdev) of each channel, the mean and standard deviation of the signal characteristic freq_COPover a two-minute moving window (20 six-second epochs). The fifth column of table311 presents current values of the freq_COP(labeled X) for the epoch immediately after the current two-minute moving window for which the statistics have been calculated. This freq_COPcharacteristic has units of Hertz (Hz) and is the frequency at the center of the bulk of signal power.
• Table311 ofFIG. 17B also includes an embodiment ofmethod elements314 and314s. After computing signal characteristics for the signals inmethod element312,method element314 sorts the channels based on the results ofmethod element312 to provide an automatic determination of which channel or channels caused the loss of timing stability detected inmethod element310. A channel-sorting standard314sis applied in the sorting ofmethod element314. Since it is assumed that the FFT-based freq_COPvalues for each channel have Gaussian distributions, one example of a useful channel-sorting standard314sis to compute Z-scores for the freq_COPvalues X for the epoch immediately after the current two-minute moving window for which the statistics have been calculated. A Z-score exceeding 2 indicates that a new freq_COP(in the set of values labeled X) has significantly (with 95% confidence) deviated from the expected distribution, thereby indicating that some change has occurred.
• The sixth, seventh and eighth columns of table311 contain absolute values of three different computed Z-scores. The isolated Z-score is computed for each channel based only on the statistics of that channel, as follows:
• 
  IsolatedZ-score=abs[(X−mean)/stdev]
• The Group Z-score is the average of the isolated Z-scores of all the channels that share a common wire index, and the combined Z-scores are the averages of the isolated and group Z-scores for each channel.
• In table311, the channel 2-HIS and channel 6-A78 signals (see grayed entries) both indicate a significant change in the current measurement of the FFT-based freq_COP; each has shifted by more than two standard deviations. (The 2-HIS isolated Z-score is 2.4, and the 6-A78 Z-score is 2.1.) However, the 6-A78 channel is on a wire that constrains the 6-A78 electrode from moving much differently than the electrodes capturing three other signals. Because the three other signals have much lower Z-scores, there is a strong indication that the channel 6-A78 isolated Z-score, especially being only a little over threshold of 2, may simply be due to natural variation and not of practical significance. (Even with a criterion like abs(Z)>2, one in20 measurements are expected to be over threshold simply by random chance.) The additional information about companion signals on the same wire (the four “A” channels) therefore assists in the recognition of a disrupted channel if the characteristic for that channel (6-A78) is near the Z>2 criterion, if many of the other companion signals less ambiguously satisfy the criterion, i.e., Z is much less than 2. Thus, in this particular example, the conclusion is that channel 2 (2-HIS) is the disrupted reference channel.
• The specific computations illustrated in tables309 and311 are not intended to be limiting with respect to the ways in which the loss of timing stability and the detection of which channel or channels have caused such loss; other specific computational approaches may be used within the scope of the inventive automatic method of determining local activation time disclosed herein.
• As indicated by the highlighted portion in gray shading in table311, cardiac electrogram channel 2-HIS is the channel whichembodiment300 of the inventive automatic method applied to this reduced channel-set example has identified as having caused the loss of timing stability. In this case, if channel 2-HIS is the current base reference channel, it will be replaced with a new reference channel selected from the remaining set of candidate reference channels, and such replacement channel is selected based on signal quality computations inmethod element316 which providemethod element318 measures of signal qualities to be compared based on a signal-quality standard318s.
• As indicated above, the qualities that may recommend a signal for use as a reference channel are: 1) high signal amplitude, 2) low signal noise, 3) low signal amplitude variability, 4) low cycle-length variability, and 5) shorter cycle length. Since several signal-quality measures determined inmethod element316 and assessed inmethod element318 using corresponding signal-quality standards318shave been described above in this document, no further examples are presented here. For the purpose selecting a new best replacement reference channel, these five channel qualities may be assessed over a longer time period since inembodiment300 of the inventive method, such assessments occur in the background, not causing any delay in the LAT mapping function. In this case, the assessment could be over a longer period such as two minutes. The most recent epoch may be incorporated into the assessment since only non-disrupted channels will be considered.
• FIG. 18 is schematic block diagram illustrating an embodiment312(2) of a first alternative method for computation of a signal characteristic for use within the method embodiment ofFIG. 15. This alternative signal characteristic computation312(2) generates a signal characteristic called epoch activity duration. Inmethod element340, values for the velocities of signal epoch data samples are computed using a method such as described above with respect toFIGS. 2-3B. The use of such method for computing signal velocities is not intended to be limiting. A number of approaches for signal velocity computation are well-known by those skilled in the art of signal processing, and such approaches are applicable within the inventive method claimed herein.
• Inmethod element342, the absolute values of the velocities are determined, and inmethod element344, an activity threshold act_THis computed. Since all MCCE signals contain some noise, it is necessary to define signal activity as occurring with signal levels above a threshold in order to avoid such noise corrupting the determination of activity duration. One useful definition of threshold act_THis four times the median of the data across the entire six-second epoch. Such threshold definition is not intended to be limiting to the present invention; other useful definitions of act_THmay be used.
• A signal is said to be active when the absolute value of the velocity is greater than act_TH. This comparison with act_THoccurs withinmethod element346 which also counts the number of signal data values which exceed act_TH.Method element348 counts the total number of signal data values, and inmethod element350, the epoch activity duration is computed as the fraction of the total signal data values exceeding act_TH.
• FIG. 18A is a plot illustrating the application of an absolute-value velocity filter to the exemplary six-second epoch cardiac channel electrogram signal ofFIG. 17A. The median of the epoch signal data ofFIG. 17A is 15 μvolts resulting in a value of act_THof 60 μvolts. The plot ofFIG. 18A indicates the value of act_THused, andFIG. 18B is a plot illustrating computation of the activity duration signal characteristic for the absolute-value velocity epoch signal ofFIG. 18A. Of the 6,000 total signal data values, 1,103 values were found to exceed act_TH=60, resulting in a value of activity duration of 0.184 for the exemplary signal epoch data ofFIG. 17A.
• It has been found that such an activity duration measurement is a useful signal characteristic by which to determine which cardiac channel or channels has/have been disturbed. Another such useful signal characteristic is simply a peak-to-peak measure across an epoch of cardiac channel data, and such a determination is described inFIG. 19 in method embodiment312(3).Method elements352 and354 find the minimum and maximum signal values, respectively, across an epoch of signal data.Method element356 computes the difference between the maximum and minimum values to yield the peak-to-peak signal characteristic.
• FIG. 19A is a table358 illustrating second alternative embodiment312(3) ofFIG. 19, showing the peak-to-peak determination of the exemplary six-second epoch of the exemplary cardiac channel electrogram signal ofFIG. 17A. The peak-to-peak value for the epoch of signal data values ofFIG. 17A is found to be 1,731 μvolts. The epoch of signal data inFIG. 17A provides a good illustration of why an epoch duration of about six seconds is desirable, and this can be seen from at least three perspectives. Six seconds is long enough to span several complete cycles of cardiac contraction; six seconds is long enough to span at least one respiratory cycle; and six seconds is long enough to obviate the need to carefully detect individual cardiac beats. By including these sources of variability entirely within, and therefore common to, an epoch, stability can be recognized from epoch to epoch.FIG. 17A clearly illustrates beat-to-beat variability within the span of a six-second epoch.
• FIG. 20 is schematic block diagram illustrating an embodiment313(4) of a third alternative method for computation of a signal characteristic for use within the method embodiment ofFIG. 15. This alternative signal characteristic computation generates a Haar-transform-based parameter called the epoch center-of-power frequency freq_COP. InFIG. 20, signal epoch data is divided into sequential segments inmethod element360. The division of signal epoch data into segments is useful for computational efficiency depending on the number of data values in a epoch of data. A Haar transform requires that the number of samples be a power of 2. Thus, for example, a six-second epoch of data captured at 1,000 sps contains 6,000 values, and such an epoch can be divided into three segments of 2,048 (2¹¹) values with only very modest overlap between the segments.
• FIGS. 20A-20C illustrate the segmenting of the exemplary data ofFIG. 17A into three sequential segments, each consisting of 2,048 values. Note that each segment is shown having sample indices from 1 to 2,048 although the actual index values range as follows: the plot ofFIG. 20A includesepoch samples 1 through 2,048; the plot ofFIG. 20B includes epoch samples 1,977 through 4,024; and the plot ofFIG. 20C includes epoch samples 3,953 through 6,000. Again, as inFIGS. 17C-17G, the abscissa of each of the three plots inFIGS. 20A-20C span from 1 through 2,048 since the Haar transform to be performed is on each segment individually although the signal data is taken from different segments of the epoch as indicated above.
• Inmethod element362, the Haar transformation coefficients for each segment of data are computed. The details of a Haar transformation are well known to those skilled in the art of signal processing. However, some details of such a computation will be described to illustrate certain aspects of method embodiment312(4).FIGS. 20D-20F are three plots of the Haar transformation coefficients for the three data segments ofFIGS. 20A-20C, respectively, resulting from the Haar transformation of such three segments of data. These coefficients are computed as illustrated in the table ofFIG. 20G. A Haar transformation is based a series of signal differences as shown inFIG. 20G, and a signal having 2,048 data values results in 2,048 Haar transformation difference-related coefficients H_i.
• The abscissa of the plots ofFIGS. 20D-20F represent the 2,048 Haar coefficients, and the plots are not frequency spectra. Unlike an FFT, the Haar transformation coefficients do not constitute a frequency spectrum but rather are a set of difference terms as described further inFIG. 20G. H₁is the sum of all 2,048 values in the segment of epoch data, and as such it represents the DC term in the signal (no differencing). H₂is the difference between the sum of the second half of the signal values minus the sum of the first half of the values, and as such, represents information in the signal related to a cycle which is the full length of the epoch segment (2.048 seconds—0.488 Hz—in the exemplary data ofFIG. 17A). The differencing continues by powers of two as shown inFIG. 20G up to H₂₀₄₈which is the difference between the last two signal values. The last 1,024 coefficients (H₁₀₂₅through H₂₀₄₈) are the differences between neighboring pairs of signal values. These coefficients relate to the highest frequency information available in the transform which in the example being discussed is 500 Hz.
• Note that in the description ofFIG. 20G, the sign convention may differ from some descriptions of the Haar transform, but since the signal characteristic being computed is related to power, the absolute value of the coefficients is used and the sign convention is not relevant. Thus, inmethod element364, an absolute value filter is applied to the Haar coefficients.
• FIG. 20H shows a table detailing the computation of a set of eleven frequency-selective aggregate magnitudes A_ifrom the 2,048 Haar transformation coefficients H_iofFIG. 20G. The Haar transformation coefficients are aggregated inmethod element366 of method embodiment312(4) ofFIG. 20 by combining the absolute values of the Haar transformation coefficients H_ihaving like time scales as shown inFIG. 20H. Thus, for example, every difference of same-time-scale sets of signal values are summed inmethod element368 of method embodiment312(4) ofFIG. 20 to form the frequency-selective aggregate magnitude related to the corresponding time scale. (The aggregate magnitudes are frequency-selective because they are formed from differences having the same time scales.)
• Such aggregation results in a spectrum-like plot of frequency-selective aggregate magnitudes A_iwhich then relate to certain frequencies which are assigned to the frequency-selective aggregate magnitudes inmethod element370 of method embodiment312(1). In the example being discussed, the frequency F_irelated to A_iis F_i=0.488×2^(i-1), and these frequencies are shown in the table ofFIG. 20H.
• The final step in method embodiment312(4), inmethod element372, is a determination of a signal characteristic based on the Haar transformation, a center-of-power frequency freq_COP. One way to determine freq_COPis to compute a first moment similar to that used for the FFT-based signal characteristic.
• 
  freq_COP=[Σ(A_i²·F_i)]/ΣA_i²
• where the sum is over the 11 frequency-selective aggregate magnitudes A_iand their corresponding frequencies F_i.
• FIGS. 20I-20K are three plots of the absolute values of the Haar transformation coefficients shown inFIGS. 20D-20F for the three data segments ofFIGS. 20A-20C, respectively, andFIGS. 20L-20N are three bar charts showing the values frequency-selective aggregate magnitudes A_ifor each of the three segments as shown inFIGS. 20A-20C and determined by applying the computations described inFIGS. 20G and 20H.
• FIG. 20P is a bar chart presenting the sum of the three sets of segments magnitudes A_ifrom the Haar transformations of the three data segments. The eleven frequency-selective aggregate magnitudes of the abscissa are used in the determination of frequency-related signal characteristic freq_COPusing the Haar-transform-based method ofFIG. 20. For the exemplary signal data ofFIG. 17A, freq_COPis found to be 141 Hz, and this is indicated in the chart ofFIG. 20P.
• Note that this frequency of 141 Hz differs from the 107 Hz found for the FFT-base computation of freq_COP. Such a difference is expected due to the logarithmic nature of the abscissa in the chart ofFIG. 20P, but since comparisons are made using the same signal characteristic computation, what is important is the ability of the signal characteristic to be sensitive to differences in the signals. All of the methods exampled above formethod element312 ofmethod embodiment300 are able to detect such differences in signals.
• While the principles of this invention have been described in connection with specific embodiments, it should be understood clearly that these descriptions are made only by way of example and are not intended to limit the scope of the invention.

Claims (27)

1. An automatic method of determining local activation time (LAT) from at least four multi-channel cardiac electrogram signals including a ventricular channel, a mapping channel and a plurality of reference channels, the method comprising:

storing the cardiac channel signals;

using the ventricular and mapping channel signals and a first reference channel signal to compute LAT values at a plurality of mapping-channel locations;

monitoring the timing stability of the first reference channel signal; and

if the timing stability of the monitored signal falls below a timing-stability standard, using a second reference channel signal to determine LAT values and avoid substantial loss of LAT values in spite of loss of timing stability.

2. The automatic LAT-determining method ofclaim 1 further including computing one or more timing offsets using pairs of the plurality of reference channel signals, a timing offset being LAT_K(J), the local activation time of a reference channel J based on a reference channel K and used to transform an LAT value based on reference channel J to an LAT value based on reference channel K.

3. The automatic LAT-determining method ofclaim 2 wherein using the second reference channel signal to determine LAT values includes transforming future LAT values such that they are based on the first reference channel.

4. The automatic LAT-determining method ofclaim 3 wherein:

LAT₂(M) is a future LAT value of mapping channel M based on the second reference channel; and

future transformed values LAT₁(M) of mapping channel M based on the first reference channel are equal to a timing offset LAT₁(2) plus LAT₂(M).

5. The automatic LAT-determining method ofclaim 2 wherein using the second reference channel signal to determine LAT values includes transforming past LAT values such that they are based on the second reference channel.

6. The automatic LAT-determining method ofclaim 5 wherein:

LAT₁(M) is a past LAT value of mapping channel M based on the first reference channel; and

past transformed values LAT₂(M) of mapping channel M based on the second reference channel are equal to a timing offset LAT₂(1) plus LAT₁(M).

7. The automatic LAT-determining method ofclaim 2 wherein the one or more timing offsets are computed at a plurality of times, and the value of each timing offset is replaced with its average over the plurality of times.

8. The automatic LAT-determining method ofclaim 7 wherein the average is computed over a predetermined number of times.

9. The automatic LAT-determining method ofclaim 2 wherein monitoring the timing stability of the first reference channel signal includes monitoring multiple timing offsets LAT₁(X) where X represents the channels with which timing offsets with the first reference channel are computed.

10. The automatic LAT-determining method ofclaim 9 further including computing a signal characteristic for the plurality of reference channels and determining therefrom which one or more channels among these reference channels has/have not lost timing stability.

11. The automatic LAT-determining method ofclaim 10 further including selecting the second reference channel signal from channels which have not lost timing stability.

12. The automatic LAT-determining method ofclaim 11 wherein selecting the second reference channel signal from channels which have not lost timing stability includes computing signal quality.

13. The automatic LAT-determining method ofclaim 10 wherein computing a signal characteristic of a signal includes computing the frequency content of the signal.

14. The automatic LAT-determining method ofclaim 13 wherein computing the frequency content of the signal includes computing a Fourier transform for a predetermined time period of the signal.

15. The automatic LAT-determining method ofclaim 13 wherein computing the frequency content of the signal includes computing a fast Fourier transform for a predetermined time period of the signal.

16. The automatic LAT-determining method ofclaim 15 wherein the computed signal characteristic is the first moment of the signal power determined from the fast Fourier transform.

17. The automatic LAT-determining method ofclaim 15 wherein the signal is segmented into a plurality of time-overlapping segment signals.

18. The automatic LAT-determining method ofclaim 17 wherein weightings are applied to each of the segment signals.

19. The automatic LAT-determining method ofclaim 18 wherein computing the fast Fourier transform of the signal includes (a) computing a signal-segment fast Fourier transform for each segment signal and (b) averaging each such signal-segment fast Fourier transform to form the fast Fourier transform of the signal.

20. The automatic LAT-determining method ofclaim 19 wherein the computed signal characteristic is the first moment of the signal power determined from the fast Fourier transform of the signal.

21. The automatic LAT-determining method ofclaim 13 wherein computing the frequency content of the signal includes computing a Haar transform for a predetermined time period of the signal.

22. The automatic LAT-determining method ofclaim 21 wherein the computed signal characteristic is the first moment of the signal power determined from the computed Haar transform.

23. The automatic LAT-determining method ofclaim 21 wherein the signal is segmented into a plurality of substantially-sequential segment signals.

24. The automatic LAT-determining method ofclaim 23 wherein computing the Haar transform of the signal includes (a) computing Haar transform coefficients for each segment signal, (b) computing absolute values of the coefficients, (c) computing a set of frequency-selective aggregate magnitudes for each segment signal by summing signal-segment Haar transform coefficients having like time scales, and (d) averaging the sets of frequency-selective aggregate magnitudes to form a single set of frequency-selective aggregate magnitudes for the signal.

25. The automatic LAT-determining method ofclaim 24 wherein the computed signal characteristic is the first moment of the signal power determined from the frequency-selective aggregate magnitudes.

26. The automatic LAT-determining method ofclaim 10 wherein the computed signal characteristic is the fraction of time within a predetermined time period of the signal at which the absolute value of signal velocity is above a predetermined threshold.

27. The automatic LAT-determining method ofclaim 10 wherein the computed signal characteristic is the maximum signal amplitude minus the minimum signal amplitude within a predetermined time period of the signal.

Images