diff --git a/lddecode/core.py b/lddecode/core.py index 7bbc1905f..1eb7592c7 100644 --- a/lddecode/core.py +++ b/lddecode/core.py @@ -2651,6 +2651,7 @@ def downscale( audio=0, final=False, lastfieldwritten=None, + shift=0 ): if lineinfo is None: lineinfo = self.linelocs @@ -2717,7 +2718,8 @@ def downscale( self.rf.downscale_sinc_lut, self.lineoffset, outwidth, - wow_level_adjust_smoothing=self.wow_level_adjust_smoothing + wow_level_adjust_smoothing=self.wow_level_adjust_smoothing, + shift=shift ) if self.rf.decode_digital_audio: diff --git a/lddecode/utils.py b/lddecode/utils.py index 0c032b437..a24ddaf2c 100644 --- a/lddecode/utils.py +++ b/lddecode/utils.py @@ -143,7 +143,7 @@ def scale(buf, begin, end, tgtlen, mult=1): @njit(nogil=True, fastmath=True) -def scale_field(buf, dsout, interpolated_pixel_locs, wowfactors, sinc_lut, lineoffset, outwidth, wow_level_adjust_smoothing = 0, level_adjust_threshold = 15): +def scale_field(buf, dsout, interpolated_pixel_locs, wowfactors, sinc_lut, lineoffset, outwidth, wow_level_adjust_smoothing = 0, level_adjust_threshold = 15, shift=0.0): # average out any unusual spikes in wow that happen on a per line basis # this indicates an hsync tbc error vs. being normal wow from playback speed variations # in this case for level adjusting we just want to fallback to the average wow to avoid a bright or dark line @@ -175,7 +175,8 @@ def scale_field(buf, dsout, interpolated_pixel_locs, wowfactors, sinc_lut, lineo level_adjust = level_adjusts[i] # reconstructs the waveform at the proper fractional sample position, undoing wow-induced timing variations - coord = np.float32(interpolated_pixel_locs[i]) + # Adding the positive shift pulls future (late) samples backward into alignment. + coord = np.float32(interpolated_pixel_locs[i] + shift) coord_int = int(coord) # fractional phase diff --git a/tests.py b/tests.py index d52596bcd..cf94c4dbb 100644 --- a/tests.py +++ b/tests.py @@ -91,10 +91,8 @@ def test_sync(filename, num_pulses=None, blank_approx=None, sync_approx=None): ldd.logger.info("test") # process.VHSDecode("infile", "outfile", ,inputfreq=samplerate_mhz, system="PAL", tape_format="VHS") - rf_options = {} - rf_options["level_detect_divisor"] = 2 rfdecoder = process.VHSRFDecode( - inputfreq=samplerate_mhz, system="PAL", tape_format="VHS", rf_options=rf_options + inputfreq=samplerate_mhz, system="PAL", tape_format="VHS" ) demod_05_data = np.loadtxt(filename) @@ -107,12 +105,13 @@ def test_sync(filename, num_pulses=None, blank_approx=None, sync_approx=None): field = FieldPALVHS(rfdecoder, data_stub) - pulses = rfdecoder.resync.get_pulses(field) + pulses = field.get_pulses() if num_pulses: assert len(pulses) == num_pulses - measured_sync, measured_blank = rfdecoder.resync._field_state.pull_levels() + measured_sync = field.sync_tip_level + measured_blank = field.blanking_level if blank_approx: assert math.isclose(measured_blank, blank_approx) @@ -122,20 +121,6 @@ def test_sync(filename, num_pulses=None, blank_approx=None, sync_approx=None): return True -def test_find_pulses(filename, num_pulses): - from vhsdecode.addons.resync import _findpulses_numba_raw - - demod_05_data = np.loadtxt(filename) - - # Just using some pre-tested values for now for optimizing function, many need changes later. - starts, lengths = _findpulses_numba_raw(demod_05_data, 3954307.8, 11.625, 1588.125) - - assert len(starts) == num_pulses - assert len(lengths) == num_pulses - assert starts[200] == 495955 - assert lengths[200] == 177 - - class SyncTest(unittest.TestCase): def test_sync_pal_good(self): blank = 4130000 @@ -151,9 +136,6 @@ def test_sync_pal_noisy(self): print("pal noisy") test_sync("PAL_NOISY.txt.gz", blank_approx=blank, sync_approx=sync) - def test_find_pulses(self): - test_find_pulses("PAL_GOOD.txt.gz", 458) - class ZCTest(unittest.TestCase): def test_calczc(self): @@ -227,10 +209,8 @@ def test_1(self): ldd.logger.info("test") demod_05_data = np.loadtxt("PAL_GOOD.txt.gz") - rf_options = {} - rf_options["level_detect_divisor"] = 2 rfdecoder = process.VHSRFDecode( - inputfreq=40, system="PAL", tape_format="VHS", rf_options=rf_options + inputfreq=40, system="PAL", tape_format="VHS" ) data_stub = {} @@ -241,11 +221,11 @@ def test_1(self): field = FieldPALVHS(rfdecoder, data_stub) - _ = rfdecoder.resync.get_pulses(field) - blank_level = rfdecoder.resync._field_state._blanklevels.current() - sync_level = rfdecoder.resync._field_state._synclevels.current() - print("blank level: ", blank_level) + _ = field.get_pulses(True) + sync_level = field.sync_tip_level + blank_level = field.blanking_level print("sync level: ", sync_level) + print("blank level: ", blank_level) if __name__ == "__main__": diff --git a/tests/run_integration.sh b/tests/run_integration.sh index 3edd202bb..38e664666 100755 --- a/tests/run_integration.sh +++ b/tests/run_integration.sh @@ -164,7 +164,6 @@ mutation _export_raw_tbc "vhs_pal" --export_raw_tbc mutation _fallback_vsync_ntsc "svhs_et_ntsc" --fallback_vsync mutation _fm_audio_notch "vhs_ntsc" --fm_audio_notch 10 mutation _ire0_adjust "vhs_ntsc" --ire0_adjust -mutation _level_detect_div5 "svhs_et_ntsc" --level_detect_divisor 5 mutation _nld "vhs_ntsc" --nld mutation _notch_ntsc "vhs_ntsc" --notch 3.58 mutation _notch_pal "vhs_pal" --notch 4.43 diff --git a/tests/unit/test_sync.py b/tests/unit/test_sync.py index 92698c660..4756c2c4e 100644 --- a/tests/unit/test_sync.py +++ b/tests/unit/test_sync.py @@ -7,16 +7,13 @@ import lddecode.core as ldd import vhsdecode.process as process from vhsdecode.field import FieldPALVHS -from vhsdecode.addons.resync import _findpulses_numba_raw - @pytest.fixture def pal_rfdecoder(): ldd.logger = logging.getLogger("test") ldd.logger.setLevel(5) - rf_options = {"level_detect_divisor": 2} return process.VHSRFDecode( - inputfreq=40, system="PAL", tape_format="VHS", rf_options=rf_options + inputfreq=40, system="PAL", tape_format="VHS" ) @@ -35,35 +32,24 @@ def _make_field(rfdecoder, filename): class TestSyncPAL: def test_sync_pal_good(self, pal_rfdecoder, data_dir): field = _make_field(pal_rfdecoder, data_dir / "PAL_GOOD.txt.gz") - pulses = pal_rfdecoder.resync.get_pulses(field) + pulses = field.get_pulses() assert len(pulses) == 458 - measured_sync, measured_blank = pal_rfdecoder.resync._field_state.pull_levels() + measured_sync, measured_blank = field.sync_tip_level, field.blanking_level assert math.isclose(measured_blank, 4133579.15, rel_tol=1e-3) assert math.isclose(measured_sync, 3840000, rel_tol=1e-3) def test_sync_pal_noisy(self, pal_rfdecoder, data_dir): field = _make_field(pal_rfdecoder, data_dir / "PAL_NOISY.txt.gz") - pal_rfdecoder.resync.get_pulses(field) - measured_sync, measured_blank = pal_rfdecoder.resync._field_state.pull_levels() - assert math.isclose(measured_blank, 4130360.76, rel_tol=1e-3) - assert math.isclose(measured_sync, 3800000, rel_tol=1e-3) - - -class TestFindPulses: - def test_find_pulses_pal_good(self, data_dir): - demod_05_data = np.loadtxt(data_dir / "PAL_GOOD.txt.gz") - starts, lengths = _findpulses_numba_raw(demod_05_data, 3954307.8, 11.625, 1588.125) - assert len(starts) == 458 - assert len(lengths) == 458 - assert starts[200] == 495955 - assert lengths[200] == 177 + _ = field.get_pulses() + measured_sync, measured_blank = field.sync_tip_level, field.blanking_level + assert math.isclose(measured_blank, 4137794.68, rel_tol=1e-3) + assert math.isclose(measured_sync, 3795981.06, rel_tol=1e-3) class TestLevelDetect: def test_level_detect_pal_good(self, pal_rfdecoder, data_dir): field = _make_field(pal_rfdecoder, data_dir / "PAL_GOOD.txt.gz") - pal_rfdecoder.resync.get_pulses(field) - blank_level = pal_rfdecoder.resync._field_state._blanklevels.current() - sync_level = pal_rfdecoder.resync._field_state._synclevels.current() - assert blank_level is not None + _ = field.get_pulses(True) + sync_level, blank_level = field.sync_tip_level, field.blanking_level assert sync_level is not None + assert blank_level is not None diff --git a/vhsdecode/addons/chromaAFC.py b/vhsdecode/addons/chromaAFC.py index 99255f80e..1fa15675a 100644 --- a/vhsdecode/addons/chromaAFC.py +++ b/vhsdecode/addons/chromaAFC.py @@ -5,8 +5,8 @@ from numpy.fft import rfft, rfftfreq import lddecode.core as ldd from scipy.signal import argrelextrema -from vhsdecode.linear_filter import FiltersClass from vhsdecode.rust_utils import sosfiltfilt_rust +from vhsdecode.linear_filter import FiltersClass twopi = 2 * np.pi diff --git a/vhsdecode/addons/resync.py b/vhsdecode/addons/resync.py deleted file mode 100644 index 55083036c..000000000 --- a/vhsdecode/addons/resync.py +++ /dev/null @@ -1,787 +0,0 @@ -from collections import namedtuple -from vhsdecode.addons.vsyncserration import VsyncSerration -import numpy as np -import vhsdecode.utils as utils -import itertools -from lddecode.utils import inrange -import lddecode.core as ldd -import math -import hashlib -from numba import njit - -# from vhsd_rust import fallback_vsync_loc_means as fallback_vsync_loc_means_r - - -@njit(cache=True, nogil=True) -def iretohz(sysparams, ire): - return sysparams.ire0 + (sysparams.hz_ire * ire) - - -@njit(cache=True, nogil=True) -def hztoire(sysparams, hz, ire0=None): - if not ire0: - ire0 = sysparams.ire0 - return (hz - ire0) / sysparams.hz_ire - - -@njit(cache=True, nogil=True) -def check_levels(data, old_sync, new_sync, new_blank, vsync_hz_ref, hz_ire, full=True): - """Check if adjusted levels are somewhat sane.""" - # ldd.logger.info("am below new blank %s , amount below half_sync %s", amount_below, amount_below_half_sync) - # ldd.logger.info("change %s _ %s", old_sync - new_sync, (hz_ire * 15)) - - blank_sync_ire_diff = (new_blank - new_sync) / hz_ire - - # ldd.logger.info("ire diff: %s", blank_sync_ire_diff) - - # Check if too far below format's standard sync, or the difference between sync and blank is too large - # to make sense - if (vsync_hz_ref - new_sync) > (hz_ire * 15) or blank_sync_ire_diff > 47: - return False - if new_sync - old_sync < (hz_ire * 5): - # Small change - probably ok - return True - - if full: - amount_below = len(np.argwhere(data < new_sync)) / len(data) - amount_below_half_sync = len(np.argwhere(data < new_blank)) / len(data) - - # If there is a lot of data below the detected vsync level, or almost no data below the detected - # 50% of hsync level it's likely the levels are not correct, so avoid adjusting. - if amount_below > 0.07 or amount_below_half_sync < 0.005: - return False - - return True - - -# search for black level on back porch -def _pulses_blacklevel(demod_05, freq_mhz: float, pulses, vsync_locs, synclevel): - if not vsync_locs or len(vsync_locs) == 0: - return None - - # take the eq pulses before and after vsync - - # We skip shorter pulses in case - before_first = vsync_locs[0] - after_last = vsync_locs[-1] - last_index = len(pulses) - 1 - - if len(vsync_locs) != 12: - # Skip pulses that are way to close together to be vsync to avoid assuming noise are pulses. - # TODO: use linelen and system for num pulses - while ( - before_first > 1 - and pulses[before_first].start - pulses[before_first - 1].start < 600 - ): - before_first -= 1 - - while ( - after_last < last_index - and pulses[after_last].start - pulses[after_last + 1].start < 600 - ): - after_last += 1 - - # TODO: This needs to be reworked for samples where the levels vary throughout the field - r1 = range(max(before_first - 5, 1), before_first) if before_first > 1 else range(0) - r2 = ( - range(after_last + 1, max(after_last + 6, last_index)) - if after_last < last_index - 1 - else range(0) - ) - - black_means = [] - - for i in itertools.chain(r1, r2): - if i < 0 or i >= len(pulses): - continue - - p = pulses[i] - if inrange(p.len, freq_mhz * 0.75, freq_mhz * 3): - mean_value = np.mean( - demod_05[int(p.start + (freq_mhz * 5)) : int(p.start + (freq_mhz * 20))] - ) - black_means.append(mean_value) - - return black_means - - -"""Pulse definition for findpulses_n. Needs to be outside the function to work with numba. -Make sure to refresh numba cache if modified. -""" -Pulse = namedtuple("Pulse", "start len") - -# Old impl for reference -# @njit(cache=True, nogil=True, parallel=True) -# def _findpulses_numba_raw_b(sync_ref, high, min_synclen, max_synclen): -# """Locate possible pulses by looking at areas within some range. -# Outputs arrays of starts and lengths -# """ -# mid_sync = high -# where_all_picture = np.where(sync_ref > mid_sync)[0] -# locs_len = np.diff(where_all_picture) -# # min_synclen = self.eq_pulselen * 1 / 8 -# # max_synclen = self.linelen * 5 / 8 -# is_sync = np.bitwise_and(locs_len > min_synclen, locs_len < max_synclen) -# where_all_syncs = np.where(is_sync)[0] -# pulses_lengths = locs_len[where_all_syncs] -# pulses_starts = where_all_picture[where_all_syncs] -# return pulses_starts, pulses_lengths - - -@njit(cache=True, nogil=True) -def _findpulses_numba_raw(sync_ref, high, min_synclen, max_synclen): - """Locate possible pulses by looking at areas within some range. - Outputs arrays of starts and lengths - """ - - in_pulse = sync_ref[0] <= high - - # Start/lengths lists - # It's possible this could be optimized further by using a different data structure here. - starts = [] - lengths = [] - - cur_start = 0 - - # Basic algorithm here is swapping between two states, going to the other one if we detect the - # current sample passed the threshold. - for pos, value in enumerate(sync_ref): - if in_pulse: - if value > high: - length = pos - cur_start - # If the pulse is in range, and it's not a starting one - if inrange(length, min_synclen, max_synclen) and cur_start != 0: - starts.append(cur_start) - lengths.append(length) - in_pulse = False - elif value <= high: - cur_start = pos - in_pulse = True - - # Not using a possible trailing pulse - # if in_pulse: - # # Handle trailing pulse - # length = len(sync_ref) - 1 - cur_start - # if inrange(length, min_synclen, max_synclen): - # starts.append(cur_start) - # lengths.append(length) - - return np.asarray(starts), np.asarray(lengths) - - -def _to_pulses_list(pulses_starts, pulses_lengths): - """Make list of Pulse objects from arrays of pulses starts and lengths""" - # Not using numba for this right now as it seemed to cause random segfault in tests. - # list(map(Pulse, pulses_starts, pulses_lengths)) - return [Pulse(z[0], z[1]) for z in zip(pulses_starts, pulses_lengths)] - - -def _findpulses_numba(sync_ref, high, min_synclen, max_synclen): - """Locate possible pulses by looking at areas within some range. - .outputs a list of Pulse tuples - """ - pulses_starts, pulses_lengths = _findpulses_numba_raw( - sync_ref, high, min_synclen, max_synclen - ) - return _to_pulses_list(pulses_starts, pulses_lengths) - - -def _fallback_vsync_loc_means( - demod_05, pulses, sample_freq_mhz: float, min_len: int, max_len: int -): - """Get the mean value of the video level inside pulses above a set threshold. - - Args: - demod_05 ([type]): Video data to get levels from. - pulses ([type]): List of detected pulses - sample_freq_mhz (float): Sample frequency of the data in mhz - min_len (int): only use pulses longer than this threshold. - max_len (int): don't use pulses longer than this threshold. - - - Returns: - [type]: a list of vsync locations in the list of pulses and a list of mean values - """ - - # ## Old impl - # vsync_locs = [] - # vsync_means = [] - - # for i, p in enumerate(pulses): - # if p.len > min_len and p.len < max_len: - # vsync_locs.append(i) - # vsync_means.append( - # np.mean( - # demod_05[ - # int(p.start + mean_pos_offset) : int( - # p.start + p.len - mean_pos_offset - # ) - # ] - # ) - # ) - - mean_pos_offset = sample_freq_mhz - - ir_pulses = [ - (i, ir_pulse) - for (i, ir_pulse) in enumerate(pulses) - if (ir_pulse.len < max_len and ir_pulse.len > min_len) - ] - - # Bail out if we didn't find matching pulses. - if not ir_pulses: - return [], [] - - def pulse_mean(locs_and_pulses): - loc, pulse = locs_and_pulses - return loc, np.mean( - demod_05[ - int(pulse.start + mean_pos_offset) : int( - pulse.start + pulse.len - mean_pos_offset - ) - ] - ) - - vsync_locs, vsync_means = list( - zip( - *map( - pulse_mean, - ir_pulses, - ) - ) - ) - - return vsync_locs, vsync_means - - -@njit(cache=True, nogil=True) -def findpulses_range(sp, vsync_hz, blank_hz=None): - """Calculate half way point between blank and sync based on ire and reference - uses calculated 0 IRE based on other params if blank_hz is not provided. - """ - if not blank_hz: - # Fall back to assume blank is at standard 0 ire. - blank_hz = iretohz(sp, 0) - sync_ire = hztoire(sp, vsync_hz) - pulse_hz_min = iretohz(sp, sync_ire - 10) - # Look for pulses at the halfway between vsync tip and blanking. - pulse_hz_max = (iretohz(sp, sync_ire) + blank_hz) / 2 - return pulse_hz_min, pulse_hz_max - - -# stores the last valid blacklevel, synclevel and vsynclocs state -# preliminary solution to fix spurious decoding halts (numpy error case) -class FieldState: - def __init__(self, sysparams): - fv = sysparams["FPS"] * 2 - ma_depth = round(fv / 5) if fv < 60 else round(fv / 6) - # ma_min_watermark = int(ma_depth / 2) - # TODO: Set to 0 for now to start using detected levels on the first field - # May want to alter later to do this more dynamically. - ma_min_watermark = 0 - self._blanklevels = utils.StackableMA( - window_average=ma_depth, min_watermark=ma_min_watermark - ) - self._synclevels = utils.StackableMA( - window_average=ma_depth, min_watermark=ma_min_watermark - ) - self._locs = None - - def set_sync_level(self, level): - self._synclevels.push(level) - - def set_levels(self, sync, blank): - self._blanklevels.push(blank) - self.set_sync_level(sync) - - def pull_sync_level(self): - return self._synclevels.pull() - - def pull_levels(self): - blevels = self._blanklevels.pull() - if blevels is not None: - return self.pull_sync_level(), blevels - else: - return None, None - - def set_locs(self, locs): - self._locs = locs - - def get_locs(self): - return self._locs - - def has_levels(self): - return self._blanklevels.has_values() and self._synclevels.has_values() - - -class Resync: - def __init__( - self, - fs, - sysparams, - sysparams_const, - thread_pool_executor, - divisor=1, - debug=False, - ): - self.divisor = divisor - self.debug = debug - self.samp_rate = fs - - if debug: - self.SysParams = sysparams.copy() - self._vsync_serration = VsyncSerration( - fs, sysparams, thread_pool_executor, divisor - ) - self._field_state = FieldState(sysparams) - self.eq_pulselen = self._vsync_serration.getEQpulselen() - self.linelen = self._vsync_serration.get_line_len() - self.use_serration = True - # This should be enough to cover all "long" pulses, - # longest variant being ones where all of vsync is just one long pulse. - self._long_pulse_max = self.linelen * 5 - # Last half-way point between blank/sync we used when looking for pulses. - self._last_pulse_threshold = findpulses_range( - sysparams_const, sysparams_const.vsync_hz - ) - - # self._temp_c = 0 - - @property - def long_pulse_max(self): - return self._long_pulse_max - - @property - def last_pulse_threshold(self): - return self._last_pulse_threshold - - def has_levels(self): - return self._field_state.has_levels() - - def _debug_field(self, sync_reference): - ldd.logger.debug( - "Hashed field sync reference %s" - % hashlib.md5(sync_reference.tobytes("C")).hexdigest() - ) - - def _pulses_blacklevel(self, field, pulses, vsync_locs, synclevel): - return _pulses_blacklevel( - field.data["video"]["demod_05"], - field.rf.freq, - pulses, - vsync_locs, - synclevel, - ) - - # checks for SysParams consistency - def _sysparams_consistency_checks(self, field): - reclamp_ire0 = False - # AGC is allowed to change two sysparams - if field.rf.useAGC: - if not self.debug: - ldd.logger.warning("Not doing consistency checks with debug disabled!") - return False - if field.rf.SysParams["ire0"] != self.SysParams["ire0"]: - ldd.logger.debug( - "AGC changed SysParams[ire0]: %.02f Hz", field.rf.SysParams["ire0"] - ) - self.SysParams["ire0"] = field.rf.SysParams["ire0"].copy() - reclamp_ire0 = True - - if field.rf.SysParams["hz_ire"] != self.SysParams["hz_ire"]: - ldd.logger.debug( - "AGC changed SysParams[hz_ire]: %.02f Hz", - field.rf.SysParams["hz_ire"], - ) - self.SysParams["hz_ire"] = field.rf.SysParams["hz_ire"].copy() - - # if self.SysParams != field.rf.SysParams: - # ldd.logger.error("SysParams changed during runtime!") - # ldd.logger.debug("Original: %s" % self.SysParams) - # ldd.logger.debug("Altered : %s" % field.rf.SysParams) - # assert False, "SysParams changed during runtime!" - - return reclamp_ire0 - - # search for sync and blanking levels from back porch - def pulses_levels( - self, field, sp, pulses, pulse_level=0, store_in_field_state=False - ): - vsync_len_px = field.usectoinpx(sp.vsync_pulse_us) - min_len = vsync_len_px * 0.8 - max_len = vsync_len_px * 1.2 - - vsync_locs, vsync_means = _fallback_vsync_loc_means( - field.data["video"]["demod_05"], pulses, field.rf.freq, min_len, max_len - ) - - # vsync_locs_2, vsync_means_2 = fallback_vsync_loc_means_r( - # field.data["video"]["demod_05"], pulses, field.rf.freq, min_len, max_len - # ) - - # print(vsync_means) - # print(vsync_means_2) - # assert np.isclose(vsync_means, vsync_means_2).all() - - if len(vsync_means) == 0: - synclevel = self._field_state.pull_sync_level() - if synclevel is None: - return None, None - else: - synclevel = np.median(vsync_means) - self._field_state.set_sync_level(synclevel) - self._field_state.set_locs(vsync_locs) - - # TODO: Think this was a bug - need to use absolute locs here, - # not position in pulse list - # if vsync_locs is None or not len(vsync_locs): - # vsync_locs = self._field_state.get_locs() - - # Now compute black level and try again - black_means = self._pulses_blacklevel(field, pulses, vsync_locs, synclevel) - # print("black_means", black_means) - - # Set to nan if empty to avoid warning. - blacklevel = ( - math.nan - if not black_means or len(black_means) == 0 - else np.median(black_means) - ) - # If black level is below sync level, something has gone very wrong. - if blacklevel < synclevel: - blacklevel = math.nan - - if False: - # store_in_field_state: - import matplotlib.pyplot as plt - - fig, (ax1, ax2) = plt.subplots(2, 1, sharex=True) - fig.set_size_inches(40, 10) - plt.text(1, 1, "black: %s sync: %s" % (blacklevel, synclevel)) - ax1.plot(field.data["video"]["demod_05"]) - ax1.axhline(synclevel) - ax1.axhline(pulse_level, color="#00FF00") - if blacklevel is not math.nan: - ax1.axhline(blacklevel, color="#000000") - for ps in pulses: - ax2.axvline(ps.start, color="#00FF00") - ax2.axvline(ps.start + ps.len, color="#0000FF") - for loc in vsync_locs: - ax2.axvline(pulses[loc].start, color="#FF0000") - # ax2.plot(self.Filters["FVideo05"]) - plt.show() - - if np.isnan(blacklevel).any() or np.isnan(synclevel).any(): - ldd.logger.debug("blacklevel or synclevel had a NaN!") - # utils.plot_scope(field.data["video"]["demod_05"], title='Failed field demod05') - sl, bl = self._field_state.pull_levels() - if bl is not None and sl is not None: - blacklevel, synclevel = bl, sl - else: - return None, None - else: - # Make sure these levels are sane before using them. - # Also don't save if we only found 1 or 2 vsyncs in case - # they were false positives. - if ( - self.level_check( - field.rf.sysparams_const, - synclevel, - blacklevel, - field.data["video"]["demod_05"], - True, - ) - and len(vsync_means) > 3 - ): - if store_in_field_state: - self._field_state.set_levels(synclevel, blacklevel) - else: - ldd.logger.debug("level check failed in pulses_levels!") - return None, None - - return synclevel, blacklevel - - def findpulses(self, sync_ref, high): - return _findpulses_numba( - sync_ref, high, self.eq_pulselen * 1 / 8, self.long_pulse_max - ) - - def _findpulses_arr(self, sync_ref, high): - return _findpulses_numba_raw( - sync_ref, high, self.eq_pulselen * 1 / 8, self.long_pulse_max - ) - - def _findpulses_arr_reduced(self, sync_ref, high, divisor, sp): - """Run findpulses using only every divisor samples""" - min_len = (self.eq_pulselen * 1 / 8) / divisor - max_len = (self.long_pulse_max) / divisor - - pulses_starts, pulses_lengths = _findpulses_numba_raw( - sync_ref[::divisor], high, min_len, max_len - ) - - pulses_starts *= divisor - pulses_lengths *= divisor - - return pulses_starts, pulses_lengths - - def add_pulselevels_to_serration_measures( - self, field, demod_05, sp, check_long=False - ): - if self._vsync_serration.has_serration(): - sync, blank = self._vsync_serration.pull_levels() - else: - # it starts finding the sync from the minima in 5 ire steps - ire_step = 5 - min_sync = np.min(demod_05) - retries = 30 - min_vsync_check = field.usectoinpx(sp.vsync_pulse_us) * 0.8 - long_pulse_min = field.usectoinpx(sp.vsync_pulse_us) * 2.6 - long_pulse_max = self.long_pulse_max - - num_assumed_vsyncs_prev = 0 - long_pulses_prev = 0 - prev_min_sync = min_sync - found_candidate = False - check_next = True - while retries > 0: - pulse_hz_min, pulse_hz_max = findpulses_range(sp, min_sync) - pulses_starts, pulses_lengths = self._findpulses_arr_reduced( - demod_05, pulse_hz_max, self.divisor, sp - ) - - # this number might need calculation - if len(pulses_lengths) > 200: - # Check that at least 2 pulses are long enough to be vsync to avoid noise - # being counted as pulses - num_assumed_vsyncs = len( - pulses_lengths[pulses_lengths > min_vsync_check] - ) - - long_pulses = 0 - - if check_long and num_assumed_vsyncs <= 2: - # If requested, we also do checks for "long" pulses found in non-standard vsync. - long_pulses = len( - pulses_lengths[ - inrange(pulses_lengths, long_pulse_min, long_pulse_max) - ] - ) - - if num_assumed_vsyncs > 4 or long_pulses >= 1: - if ( - num_assumed_vsyncs == 12 or long_pulses == 2 - ) and not check_next: - # if we have exactly 12 vsyncs meaning we likely found all of them and no more - # (or 2 long pulses) we are likely good. - break - elif ( - not found_candidate - or num_assumed_vsyncs > num_assumed_vsyncs_prev - or long_pulses > long_pulses_prev - ): - # We found a set with at least some vsyncs or long pulses... - found_candidate = True - num_assumed_vsyncs_prev = num_assumed_vsyncs - long_pulses_prev = long_pulses - prev_min_sync = min_sync - # ...so do one more iteration to see if we get closer to the expected number, - # if not use the currently found levels. - check_next = True - elif ( - num_assumed_vsyncs < num_assumed_vsyncs_prev - or long_pulses < long_pulses_prev - or check_next is False - ): - # We found less so go back to previous levels and end the search. - min_sync = prev_min_sync - pulse_hz_min, pulse_hz_max = findpulses_range(sp, min_sync) - pulses_starts, pulses_lengths = self._findpulses_arr( - demod_05, pulse_hz_max - ) - break - else: - check_next = False - - min_sync = iretohz(sp, hztoire(sp, min_sync) + ire_step) - retries -= 1 - - pulses = _to_pulses_list(pulses_starts, pulses_lengths) - - sync, blank = self.pulses_levels(field, sp, pulses, pulse_hz_max) - # chewed tape case - if sync is None or blank is None: - ldd.logger.debug("Level detection failed - sync or blank is None") - return - - # the tape chewing test passed, then it should find sync - pulse_hz_min, pulse_hz_max = findpulses_range(sp, sync, blank_hz=blank) - - pulses = self.findpulses(demod_05, pulse_hz_max) - - if False: - import matplotlib.pyplot as plt - - # ldd.logger.info("hz to ire %s", hztoire(sp, blank)) - fig, (ax1, ax2) = plt.subplots(2, 1, sharex=True) - ax1.plot(field.data["video"]["demod_05"]) - ax1.axhline(sync) - ax1.axhline(pulse_hz_max, color="#00FF00") - # ax1.axhline((sync + blank) / 2.0 , color="#FF0000") - ax1.axhline(blank, color="#000000") - for p in pulses: - ax2.axvline(p.start, color="#00FF00") - ax2.axvline(p.start + p.len, color="#0000FF") - # ax2.plot(self.Filters["FVideo05"]) - plt.show() - - f_sync, f_blank = self.pulses_levels( - field, sp, pulses, pulse_hz_max, store_in_field_state=True - ) - if f_sync is not None and f_blank is not None: - self._vsync_serration.push_levels((f_sync, f_blank)) - else: - ldd.logger.debug( - "Level detection had issues, so don't store anything in VsyncSerration." - ) - - # Do a level check - def level_check(self, sysparams_const, sync, blank, sync_reference, full=True): - vsync_hz = ( - sysparams_const.vsync_hz - ) # field.rf.iretohz(field.rf.SysParams["vsync_ire"]) - # TODO: See if we need to read vsync_ire from sysparams here - return check_levels( - sync_reference, - vsync_hz, - sync, - blank, - sysparams_const.vsync_hz, - sysparams_const.hz_ire, - full, - ) - - def get_pulses(self, field, check_levels=True): - if self.use_serration: - return self._get_pulses_serration(field, check_levels) - else: - import vhsdecode.leveldetect - - sync, blank = None, None - if self._field_state.has_levels(): - sync, blank = self._field_state.pull_levels() - pulses = self.findpulses( - field.data["video"]["demod_05"], (blank + sync) / 2 - ) - if len(pulses) > 200 and len(pulses) < 800: - return pulses - ldd.logger.info("Re-checking levels..") - - def_sync = field.rf.iretohz(field.rf.SysParams["vsync_ire"]) - def_blank = field.rf.iretohz(field.rf.SysParams["ire0"]) - sync, blank = vhsdecode.leveldetect.find_sync_levels( - field.data["video"]["demod_05"], - def_sync, - def_blank, - field.get_linefreq(), - ) - self._field_state.set_levels(sync, blank) - - return self.findpulses(field.data["video"]["demod_05"], (blank + sync) / 2) - - def _get_pulses_serration(self, field, check_levels): - """Find sync pulses in the demodulated video signal""" - - sp = field.rf.sysparams_const - - sync_reference = field.data["video"]["demod_05"] - if self.debug: - self._debug_field(sync_reference) - - if check_levels or not self._field_state.has_levels(): - # measures the serration levels if possible - # Just skip the serration stuff on 405 and 819 line stuff for now - # as it won't work. - if not (field.rf.color_system == "405" or field.rf.color_system == "819"): - self._vsync_serration.work(sync_reference) - # adds the sync and blanking levels from the back porch - self.add_pulselevels_to_serration_measures( - field, sync_reference, sp, field.rf.options.fallback_vsync - ) - - # safe clips the bottom of the sync pulses but leaves picture area unchanged - # NOTE: Disabled for now as it doesn't seem to have much purpose at the moment and can - # cause weird artifacts on the output. - demod_data = ( - field.data["video"]["demod"] - # if not field.rf.options.sync_clip - # else self._vsync_serration.safe_sync_clip( - # sync_reference, field.data["video"]["demod"] - # ) - ) - - # if self._temp_c == 1: - # np.savetxt("PAL_GOOD.txt.gz", sync_reference) - # self._temp_c += 1 - - # if it has levels, then compensate blanking bias - if self._vsync_serration.has_levels() or self._field_state.has_levels(): - if self._vsync_serration.has_levels(): - new_sync, new_blank = self._vsync_serration.pull_levels() - if self.level_check(sp, new_sync, new_blank, sync_reference): - sync, blank = new_sync, new_blank - elif self._field_state.has_levels(): - sync, blank = self._field_state.pull_levels() - ldd.logger.debug( - "Level check failed on serration measured levels [new_sync: %s, new_blank: %s], falling back to levels from FieldState [sync %s, blank %s].", - new_sync, - new_blank, - sync, - blank, - ) - else: - # Check failed on serration levels and field state does not contain levels - # TODO: Handle properly if this occurs - ldd.logger.debug( - "Level check failed on serration measured levels, using defaults." - ) - - sync = sp.ire0 - blank = sp.vsync_hz - else: - sync, blank = self._field_state.pull_levels() - - if self._sysparams_consistency_checks(field): - field.rf.SysParams["ire0"] = blank - - dc_offset = sp.ire0 - blank - sync_reference += dc_offset - if not field.rf.options.disable_dc_offset: - demod_data += dc_offset - field.data["video"]["demod"] = demod_data - sync, blank = sync + dc_offset, blank + dc_offset - - field.data["video"]["demod_05"] = sync_reference - pulse_hz_min, pulse_hz_max = findpulses_range(sp, sync) - else: - # pass one using standard levels (fallback sync logic) - # pulse_hz range: vsync_ire - 10, maximum is the 50% crossing point to sync - pulse_hz_min, pulse_hz_max = findpulses_range(sp, sp.vsync_hz) - - # checks if the DC offset is abnormal before correcting it - new_sync = self._vsync_serration.mean_bias() - vsync_hz = sp.vsync_hz - new_blank = iretohz(sp, hztoire(sp, new_sync) / 2) - - check = self.level_check(sp, new_sync, new_blank, sync_reference) - if ( - not field.rf.options.disable_dc_offset - and not pulse_hz_min < new_sync < vsync_hz - and check - ): - field.data["video"]["demod_05"] = sync_reference - new_sync + vsync_hz - field.data["video"]["demod"] = demod_data - new_sync + vsync_hz - - self._last_pulse_threshold = pulse_hz_max - - return self.findpulses(field.data["video"]["demod_05"], pulse_hz_max) diff --git a/vhsdecode/addons/vsyncserration.py b/vhsdecode/addons/vsyncserration.py deleted file mode 100644 index 0a8aeffe0..000000000 --- a/vhsdecode/addons/vsyncserration.py +++ /dev/null @@ -1,433 +0,0 @@ -# It search for the VBI serration EQ pulses, -# locates them, and extracts the sync and blanking levels. -# Also gives the base for level clamping, and maybe genlocking/vroom prevention - -from vhsdecode.utils import ( - firdes_lowpass, - firdes_highpass, - plot_scope, - dualplot_scope, - zero_cross_det, - StackableMA, -) - -from vhsdecode.linear_filter import FiltersClass, chainfiltfilt_b -import numpy as np -from os import getpid -from numba import njit - -import lddecode.core as ldd - - -# from frequency to samples -def f_to_samples(samp_rate, frequency): - return samp_rate / frequency - - -# from time to samples -def t_to_samples(samp_rate, time): - return f_to_samples(samp_rate, 1 / time) - - -# applies a filtfilt to the data over the array of filters -def _chainfiltfilt(data, filters): - for filt in filters: - data = filt.filtfilt(data) - return data - - -def _local_minima_indices(data): - """Return local minima indexes using a low-allocation vectorized comparison.""" - if len(data) < 3: - return np.empty(0, dtype=np.int64) - - mid = data[1:-1] - prev = data[:-2] - nxt = data[2:] - - # Keep allocations minimal: build one boolean mask and refine it in-place. - minima_mask = np.empty(len(mid), dtype=bool) - np.less(mid, prev, out=minima_mask) - np.less_equal(mid, nxt, out=minima_mask, where=minima_mask) - - return np.nonzero(minima_mask)[0] + 1 - - -# clips any sync pulse that satisfies min_synclen < pulse_len < max_synclen -# @njit(cache=True) -# def _safe_sync_clip(sync_ref, data, levels, eq_pulselen): -# sync, blank = levels -# mid_sync = (sync + blank) / 2 -# where_all_picture = np.where(sync_ref > mid_sync)[0] -# locs_len = np.diff(where_all_picture) -# min_synclen = eq_pulselen * 3 / 4 -# max_synclen = eq_pulselen * 3 -# is_sync = np.bitwise_and(locs_len > min_synclen, locs_len < max_synclen) -# where_all_syncs = np.where(is_sync)[0] -# clip_from = where_all_picture[where_all_syncs] -# clip_len = locs_len[where_all_syncs] -# for ix, begin in enumerate(clip_from): -# data[begin : begin + clip_len[ix]] = sync -# -# return data - - -# encapsulates the serration search logic -class VsyncSerration: - def __init__( - self, - fs, - sysparams, - thread_pool_executor, - divisor=1, - show_decoded_serration=False, - ): - self._divisor = divisor - self.show_decoded = show_decoded_serration - self.samp_rate = fs / self._divisor - self._thread_pool_executor = thread_pool_executor - fv = sysparams["FPS"] * 2 - fh = sysparams["FPS"] * sysparams["frame_lines"] - - # parameter, harmonic limit of the envelope search (with respect of vertical frequency) - venv_limit = 5 - # parameter, divisor of fh for limiting the bandwidth of power_ratio_search() - serration_limit = 3 - # parameter, depth/window of the moving averaging - ma_depth = 2 - ma_min_watermark = 1 - - # used on vsync_envelope_simple() (search for video amplitude pinch) - iir_vsync_env = firdes_lowpass(self.samp_rate, fv * venv_limit, 1e3) - self.vsyncEnvFilter = FiltersClass( - iir_vsync_env[0], iir_vsync_env[1], self.samp_rate - ) - - # used in power_ratio_search(), it makes a bandpass filter - # cannot design it as a bandpass with the given constraints - iir_serration_base_lo = firdes_highpass(self.samp_rate, fh, fh) - iir_serration_base_hi = firdes_lowpass(self.samp_rate, fh, fh) - self.serrationFilter_base = [ - FiltersClass( - iir_serration_base_lo[0], iir_serration_base_lo[1], self.samp_rate - ), - FiltersClass( - iir_serration_base_hi[0], iir_serration_base_hi[1], self.samp_rate - ), - ] - - iir_serration_envelope_lo = firdes_lowpass( - self.samp_rate, fh / serration_limit, fh / 2 - ) - - self.serrationFilter_envelope = FiltersClass( - iir_serration_envelope_lo[0], iir_serration_envelope_lo[1], self.samp_rate - ) - - # -- end of uses of power_ratio_search() - - # several timing related constants - self.eq_pulselen = round( - t_to_samples(self.samp_rate, sysparams["eqPulseUS"] * 1e-6) - ) - self.vsynclen = round(f_to_samples(self.samp_rate, fv)) - self.linelen = round(f_to_samples(self.samp_rate, fh)) - line_time = 1 / fh - vbi_time = 6.5 * line_time - self.vbi_time_range = t_to_samples( - self.samp_rate, vbi_time * 3 / 4 - ), t_to_samples(self.samp_rate, vbi_time * 5 / 4) - - # result storage instances - self.levels = StackableMA( - window_average=ma_depth, min_watermark=ma_min_watermark - ), StackableMA( - window_average=ma_depth, min_watermark=ma_min_watermark - ) # sync, blanking - - self.sync_level_bias = np.array([]) - self.fieldcount = 0 - self.pid = getpid() - self.found_serration = False - - def getEQpulselen(self): - return self.eq_pulselen * self._divisor - - def get_line_len(self): - return self.linelen * self._divisor - - # returns the measured sync level and blank level - def pull_levels(self): - sync, blank = self.levels[0].pull(), self.levels[1].pull() - return sync, blank - - # returns true if it has levels above the min_watermark - def has_levels(self): - return self.levels[0].has_values() and self.levels[1].has_values() - - def has_serration(self): - return self.found_serration - - # it adds external levels to the stack - def push_levels(self, levels): - for ix, level in enumerate(levels): - self.levels[ix].push(level) - - # only used when printing the charts - def _mutemask(self, raw_locs, blocklen, pulselen): - mask = np.zeros(blocklen) - locs = raw_locs[np.where(raw_locs < blocklen - pulselen)[0]] - for loc in locs: - mask[loc : loc + pulselen] = [1] * pulselen - return mask[:blocklen] - - # this may need tweak - def _vsync_envelope_simple(self, data): - # hi_part = np.clip(data, a_max=np.max(data), a_min=0) - hi_part = data - hi_filtered = self.vsyncEnvFilter.filtfilt(hi_part) - return hi_filtered, np.min(data) - - # does vsync_envelope_simple in forward and reverse direction, - # then assembles both halves as one result. - # It is a hack to avoid edge distortion when using lowpass filters - # of very low cutoff - def _vsync_envelope_double(self, data): - half = int(len(data) / 2) - - hi_part = np.clip(data, a_max=None, a_min=0) - - def vsync_env_rev(hp): - return np.flip(self.vsyncEnvFilter.filtfilt(np.flip(hp))) - - # Do the two filter operations on separate threads for a speedup. - forward_f = self._thread_pool_executor.submit( - self._vsync_envelope_simple, hi_part - ) - # reverse_t_f = self._thread_pool_executor.submit(self.vsync_envelope_simple, np.flip(hi_part)) - reverse_t_f = self._thread_pool_executor.submit(vsync_env_rev, hi_part) - forward = forward_f.result() - reverse = reverse_t_f.result() - # # Non-threaded version: - # b_forward = self.vsync_envelope_simple(hi_part) - # b_reverse_t = self.vsync_envelope_simple(np.flip(hi_part)) - # b_reverse = np.flip(b_reverse_t[0]), b_reverse_t[1] - - # end of forward + beginning of reverse - # Re-use existing array instead of allocating a new one. - # Maybe there's a better way to do this. - result_temp = hi_part - result_temp[:half] = reverse[:half] - result_temp[half:] = forward[0][half:] - result = result_temp, forward[1] - # # Old version - # b_result = (np.append(b_reverse[0][:half], b_forward[0][half:]), b_forward[1]) - - # dualplot_scope(forward[0], reverse[0]) - # dualplot_scope(result[0], result[1], title="VBI envelope") - return result - - # measures the harmonics of the EQ pulses - def _power_ratio_search(self, data): - try: - first_harmonic = chainfiltfilt_b(data, self.serrationFilter_base) - # Make sure we do this in place to avoid allocating an extra array. - first_harmonic **= 2 - first_harmonic = self.serrationFilter_envelope.filtfilt(first_harmonic) - return _local_minima_indices(first_harmonic) - except MemoryError: - ldd.logger.warning( - "Low-memory condition during serration harmonic minima search; using fallback path." - ) - return np.empty(0, dtype=np.int64) - - # fills in missing VBI positions when possible - @staticmethod - @njit(nogil=True, cache=True, fastmath=True) - def _vsync_arbitrage(vsynclen, where_allmin, serrations, datalen): - result = [] - if len(where_allmin) > 1: - valid_serrations = [] - - # select serration - for id, edge in enumerate(serrations): - for s_min in where_allmin: - next_serration_id = min(id + 1, len(serrations) - 1) - if edge <= s_min <= serrations[next_serration_id]: - valid_serrations.append(edge) - - for serration in valid_serrations: - if serration - vsynclen >= 0 or serration + vsynclen <= datalen - 1: - result.append(serration) - elif len(where_allmin) == 1: - if where_allmin[0] + vsynclen < datalen - 1: - result.append(where_allmin[0]) - result.append(where_allmin[0] + vsynclen) - else: - result.append(where_allmin[0]) - result.append(max(where_allmin[0] - vsynclen, 0)) - else: - result = None - - return result - - # extracts the level from a valid serration - @staticmethod - @njit(nogil=True, cache=True, fastmath=True) - def _get_serration_sync_levels(serration): - half_amp = np.mean(serration) - peaks = np.where(serration > half_amp)[0] - valleys = np.where(serration <= half_amp)[0] - levels = np.median(serration[valleys]), np.median(serration[peaks]) - return levels - - # validates the found section as a serration - def _search_eq_pulses(self, data, pos, linespan=30): - start, end = max(0, pos - self.linelen * linespan), min( - len(data) - 1, pos + self.linelen * linespan - ) - min_block = data[start:end] - min_block_min = np.min(min_block) - level = (np.median(min_block) - min_block_min) / 2 - level += min_block_min - zero_block = min_block - level - sync_pulses = zero_cross_det(zero_block) - diff_sync = np.diff(sync_pulses) - - where_min_diff = np.where( - np.logical_and( - self.eq_pulselen * 0.2 < diff_sync, diff_sync < self.eq_pulselen * 5 / 4 - ) - )[0] - - # a valid serration should count about 11 pulses, but dropouts and noise - if 9 <= len(where_min_diff) <= 12: - eq_s, eq_e = sync_pulses[where_min_diff[0]], min( - int(sync_pulses[where_min_diff[-1:][0]] + self.eq_pulselen / 2), - len(data) - 1, - ) - data_s, data_e = eq_s + start, eq_e + start - serration = data[data_s:data_e] - - # validates it by time length, (original version 17e3 and 23e3) - # now calculated at initialization - if self.vbi_time_range[0] < len(serration) < self.vbi_time_range[1]: - self.found_serration = True - self.push_levels(VsyncSerration._get_serration_sync_levels(serration)) - - if self.show_decoded: - sync, blank = self.pull_levels() - marker = np.ones(len(serration)) * blank - dualplot_scope( - serration, - marker, - title="VBI EQ serration + measured blanking level", - ) - return True, data_s, data_e - else: - return False, None, None - else: - return False, None, None - - def mean_bias(self): - return np.mean(self.sync_level_bias) - - def _remove_bias(self, data): - return data - self.sync_level_bias - - # this is the start-of-search - def _vsync_envelope(self, data, padding=1024): # 0x10000 - padded = np.append(np.flip(data[:padding]), data) - forward = self._vsync_envelope_double(padded) - self.sync_level_bias = forward[1] - - # Optimization - do in place - # diff = np.add(forward[0][padding:], -self.sync_level_bias) - forward[0][padding:] -= self.sync_level_bias - diff = forward[0][padding:] - # Since we modified this, make sure it's not re-used accidentaly. - del forward - try: - where_allmin = _local_minima_indices(diff) - except MemoryError: - ldd.logger.warning( - "Low-memory condition during video envelope minima search; using fallback logic." - ) - where_allmin = np.empty(0, dtype=np.int64) - if len(where_allmin) > 0: - serrations = self._power_ratio_search(padded) - where_min = VsyncSerration._vsync_arbitrage( - self.vsynclen, where_allmin, serrations, len(padded) - ) - serration_locs = list() - if len(where_min) > 0: - mask_len = self.linelen * 5 - state = False - tasks = [ - self._thread_pool_executor.submit( - self._search_eq_pulses, data, w_min - ) - for w_min in where_min - ] - - for task in tasks: - state, serr_loc, serr_len = task.result() - if state: - serration_locs.append(serr_loc) - mask_len = serr_len - serr_loc - - if self.show_decoded: - data_copy = data.copy() # self.remove_bias(data) - if len(serration_locs) > 0: - mask = self._mutemask( - np.array(serration_locs), len(data_copy), mask_len - ) - dualplot_scope( - data_copy, - np.clip( - mask * np.amax(data_copy), - a_max=np.amax(data_copy), - a_min=np.amin(data_copy), - ), - title="VBI position", - ) - else: - plot_scope(data_copy, title="Missing serration measure") - ldd.logger.warning("A serration measure is missing") - return state - else: - # dualplot_scope(forward[0], forward[1], title='unexpected arbitrage') - ldd.logger.warning("Unexpected vsync arbitrage") - return None - else: - # dualplot_scope(forward[0], forward[1], title='unexpected, there is no minima') - ldd.logger.warning("Unexpected video envelope") - return None - - # this runs the measures - def work(self, data): - self.found_serration = False - self._vsync_envelope(data[:: self._divisor]) - if self.has_levels() and self.found_serration: - ldd.logger.debug( - "VBI serration levels %d - Sync tip: %.02f kHz, Blanking (ire0): %.02f kHz" - % ( - self.levels[0].size(), - self.pull_levels()[0] / 1e3, - self.pull_levels()[1] / 1e3, - ) - ) - elif self.fieldcount % 10 == 0: - ldd.logger.debug( - "VBI EQ serration pulses search failed (using fallback logic)" - ) - - self.fieldcount += 1 - - # safe clips the bottom of the sync pulses, but not the picture area - # def safe_sync_clip(self, sync_ref, data): - # if self.has_levels(): - # data = _safe_sync_clip( - # sync_ref, data, self.pull_levels(), self.getEQpulselen() - # ) - # return data diff --git a/vhsdecode/chroma.py b/vhsdecode/chroma.py index 53e0b6e6f..d2919dfb8 100644 --- a/vhsdecode/chroma.py +++ b/vhsdecode/chroma.py @@ -9,6 +9,7 @@ import numba from numba import njit from numba.experimental import jitclass +from functools import cache @njit(cache=True, nogil=True, fastmath=True) @@ -27,70 +28,106 @@ def chroma_to_u16(chroma): return out @njit(cache=False, nogil=True, fastmath=True) -def chroma_automatic_gain(chroma, burst_abs_ref, phase_sequence, burst_detected_line, smoothing_window=8, k=2.0): +def chroma_automatic_gain( + chroma, + burst_abs_ref, + phase_sequence, + burst_detected_line, + sync_tip_len, + smoothing_window=8, + k=2.0 +): burst_count = len(phase_sequence) - raw_gains = np.zeros(burst_count, dtype=np.float64) - # extract gain values - valid_gains_list = [] + raw_gains = np.empty(burst_count, dtype=np.float64) + valid_gains = np.empty(burst_count, dtype=np.float64) + valid_amps = np.empty(burst_count, dtype=np.float64) + valid_count = 0 + + # extract gain values and track valid amplitudes for i in range(burst_count): current_burst = phase_sequence[i] - curr_amp = current_burst.amplitude if current_burst.amplitude != 0 else 1e-8 - raw_gains[i] = burst_abs_ref / curr_amp - - if current_burst.line_number >= burst_detected_line: - valid_gains_list.append(raw_gains[i]) + current_amp = current_burst.amplitude if current_burst.amplitude != 0 else 1e-8 + + raw_gain = burst_abs_ref / current_amp + raw_gains[i] = raw_gain - # calculate MAD - if len(valid_gains_list) > 0: - valid_gains = np.array(valid_gains_list) - median_gain = np.median(valid_gains) - mad_gain = np.median(np.abs(valid_gains - median_gain)) + if current_burst.line_number >= burst_detected_line: + valid_gains[valid_count] = raw_gain + valid_amps[valid_count] = current_amp + valid_count += 1 + + # calculate MAD threshold for gain adjustment + if valid_count > 0: + active_gains = valid_gains[:valid_count] + median_gain = np.median(active_gains) + mad_gain = np.median(np.abs(active_gains - median_gain)) max_allowable_gain = median_gain + (k * mad_gain) else: max_allowable_gain = 1.0 + # clamp gains + clamped_gains = np.minimum(raw_gains, max_allowable_gain) + # calculate smoothing - smoothed_gains = np.zeros(burst_count, dtype=np.float64) + smoothed_gains = np.empty(burst_count, dtype=np.float64) half_w = smoothing_window // 2 for i in range(burst_count): - # Apply moving average directly to the clamped values start = max(0, i - half_w) end = min(burst_count, i + half_w + 1) - - window_sum = 0.0 - for w in range(start, end): - window_sum += min(raw_gains[w], max_allowable_gain) - smoothed_gains[i] = window_sum / (end - start) + smoothed_gains[i] = np.sum(clamped_gains[start:end]) / (end - start) - # apply gain + # apply gain, calculate noise floor + noise_sum = 0 + noise_samples = 0 for i in range(burst_count): current_burst = phase_sequence[i] current_burst_start = current_burst.start - next_burst_start = phase_sequence[i + 1].start if i < burst_count - 1 else len(chroma) + + if i < burst_count - 1: + next_burst_start = phase_sequence[i + 1].start + else: + next_burst_start = len(chroma) if current_burst.line_number < burst_detected_line: - chroma[current_burst_start:next_burst_start] = 0 + chroma[current_burst_start:next_burst_start] = 0.0 else: gain_start = smoothed_gains[i] - gain_end = smoothed_gains[i + 1] if i < burst_count - 1 else smoothed_gains[i] - gain_increment = (gain_end - gain_start) / (next_burst_start - current_burst_start) - gain = gain_start + if i < burst_count - 1: + gain_end = smoothed_gains[i + 1] + else: + gain_end = smoothed_gains[i] + + length = next_burst_start - current_burst_start + if length > 0: + gain_increment = (gain_end - gain_start) / length + gain = gain_start + + # apply gain + for j in range(current_burst_start, next_burst_start): + chroma[j] = chroma[j] * gain + gain += gain_increment + + # get noise floor of sync tip area in the chroma channel + sync_tip = chroma[next_burst_start + 4 - sync_tip_len : next_burst_start - 4] + + # MAD + med = np.median(sync_tip) + abs_dev = np.abs(sync_tip - med) + mad = np.median(abs_dev) - for j in range(current_burst_start, next_burst_start): - chroma[j] = chroma[j] * gain - gain += gain_increment + # Accumulate the noise floor + noise_sum += mad * 1.4826 + noise_samples += 1 + + # Calculate the average noise floor for the entire processed region + noise_floor = noise_sum / noise_samples if noise_samples > 0 else 0.0 - # calculate RMS - if len(valid_gains_list) > 0: - sq_sum = 0.0 - for g in valid_gains_list: - # Reconstruct original amplitude from valid gains to calculate signal RMS - amp = burst_abs_ref / g - sq_sum += amp * amp - return np.sqrt(sq_sum / len(valid_gains_list)) + # return mean of means + if valid_count > 0: + return np.mean(valid_amps[:valid_count]), noise_floor - return 0.0 + return 0.0, noise_floor @njit(cache=True, nogil=True) @@ -775,7 +812,6 @@ def upconvert_chroma( @njit(nogil=True, cache=False, fastmath=False) def upconvert_chroma_phase_comp( chroma, - uphet, lineoffset, outwidth, phase_rotation_sequence, @@ -840,7 +876,6 @@ def upconvert_chroma_phase_comp( # Slice target and source arrays to provide direct contiguous memory views chroma_slice = chroma[linestart:lineend] - uphet_slice = uphet[linestart:lineend] # No outer dependencies so LLVM can vectorize for k in range(outwidth): @@ -848,7 +883,7 @@ def upconvert_chroma_phase_comp( theta_k = theta_0 + alpha * local_idx[k] + beta * (local_idx[k] * local_idx[k]) # Vectorized trigonometric and arithmetic execution - uphet_slice[k] = chroma_slice[k] * -np.cos(theta_k) - dc_val + chroma_slice[k] = chroma_slice[k] * -np.cos(theta_k) - dc_val @njit(cache=True, nogil=True) @@ -964,7 +999,7 @@ def decode_chroma_phase_rotation( burstarea, field.rf.fsc_wave, field.rf.fsc_cos_wave, - field.rf.chroma_afc.fsc_mhz * 1e6, + field.rf.SysParams['fsc_mhz'] * 1e6, detect_chroma_track_phase, rotation_check_start_line, # check for track phase rotation around the headswitching area (bottom of field) field.rf.options.enable_color_killer, @@ -1667,6 +1702,85 @@ def _process_chroma_secam_method1(field, chroma, linesout, outwidth, burstarea): return uphet +@cache +def _gen_chroma_fft_filter( + filter_len: int, + fsc: float, + color_under_carrier_f: float, + bw_lower_hz: float, + heterodyne_attenuation_db: float, + order: int, +) -> np.ndarray: + """ + Generate asymmetric Super-Gaussian bandpass mask in the frequency domain. + """ + + # calculate upper bandwidth limit that is required to remove the heterodyne up conversion product + # at the supplied attenuation + A = 10.0 ** (-abs(heterodyne_attenuation_db) / 20.0) + delta_f = 2.0 * color_under_carrier_f + exponent = 1.0 / (2.0 * order) + bw_upper_hz = delta_f / ((-np.log(A)) ** exponent) + + freqs_up = sps_fft.rfftfreq(filter_len, d=1.0 / (fsc * 4.0)) + mask = np.zeros_like(freqs_up, dtype=np.float64) + + lower_idx = freqs_up <= fsc + upper_idx = freqs_up > fsc + + # asymmetric Super-Gaussian evaluation around subcarrier (fsc) + mask[lower_idx] = np.exp(-((freqs_up[lower_idx] - fsc) / bw_lower_hz) ** (2 * order)) + mask[upper_idx] = np.exp(-((freqs_up[upper_idx] - fsc) / bw_upper_hz) ** (2 * order)) + + return mask + + +def filter_chroma_fft( + uphet: np.ndarray, + fsc: float, + color_under_carrier_f: float, + bw_lower_hz: float, # Lower chroma bandwidth (1.3 MHz below fsc) + heterodyne_attenuation_db: float, # Rejection target at sum product (dB) + order: int = 2, # filter order + pad_samples: int = 256, +) -> np.ndarray: + """ + Zero-phase FFT bandpass filter using a Super-Gaussian mask. + Dynamically computes bw_upper_hz from color_under_carrier_f to guarantee + stopband attenuation at the heterodyne sum product. + """ + N_raw = len(uphet) + + # pad to nearest fast FFT length (combines 2, 3, 5, 7 prime factors) + min_pad_len = N_raw + 2 * max(1, int(pad_samples)) + N_up = sps_fft.next_fast_len(min_pad_len) + + # Symmetric pad calculation to center the signal in the fast length array + pad_left = (N_up - N_raw) // 2 + pad_right = N_up - N_raw - pad_left + + x_padded = np.pad(uphet, (pad_left, pad_right), mode='reflect') + + mask = _gen_chroma_fft_filter( + N_up, + fsc, + color_under_carrier_f, + bw_lower_hz, + heterodyne_attenuation_db, + order + ) + + # apply filter against Forward Real FFT + F_filtered = sps_fft.rfft(x_padded) + F_filtered *= mask + + # return to real signal + chroma_padded = sps_fft.irfft(F_filtered, n=N_up) + + # remove padding + return chroma_padded[pad_left : pad_left + N_raw] + + def process_chroma( field, disable_deemph=False, @@ -1678,15 +1792,30 @@ def process_chroma( linesout = field.outlinecount outwidth = field.outlinelen - uphet = np.zeros((linesout * outwidth), dtype=np.float32) if field.burst_detected_line == -1: # skip chroma if the color killer is active for the whole field - return uphet + return np.zeros((linesout * outwidth), dtype=np.float32) + + if ( + not field.rf.options.disable_phase_correction + and field.rf.color_system == "NTSC" + ): + field.fieldPhaseID, target_phase = ntsc_color_framing_map[ + (field.isFirstField, (field.field_number // 2) % 2) + ] + chroma_shift_direction = 1 if target_phase else -1 + else: + chroma_shift_direction = 0 # Run TBC/downscale on chroma (if new field, else uses cache) # Cached if chroma process is run multiple times on one field due to track detection. if field.chroma_tbc_buffer is None: - chroma, _, _ = ldd.Field.downscale(field, channel="demod_burst") + # shift the chroma to reverse group delay caused by the color under heterodyne filter + # this is dependent on color framing, and is disabled if color framing is disabled + # TODO: shift amount may need tuning / needs validation + chroma_subcarrier_delay_cycles = field.rf.SysParams['fsc_mhz'] * 1e6 / (2.0 * np.pi * field.rf.DecoderParams["color_under_carrier"]) + chroma_subcarrier_delay_samples = chroma_subcarrier_delay_cycles * 4 + chroma, _, _ = ldd.Field.downscale(field, channel="demod_burst", shift=chroma_subcarrier_delay_samples * chroma_shift_direction) # If chroma AFC is enabled if field.rf.do_cafc: @@ -1723,7 +1852,7 @@ def process_chroma( outwidth, porch_window, field.rf.chroma_afc.true_samp_rate, - field.rf.chroma_afc.color_under, + field.rf.DecoderParams["color_under_carrier"], ) if carrier_offset is not None: field.rf.secam_servo_avg.push(carrier_offset) @@ -1755,9 +1884,6 @@ def process_chroma( not field.rf.options.disable_phase_correction and field.rf.color_system == "NTSC" ): - field.fieldPhaseID, target_phase = ntsc_color_framing_map[ - (field.isFirstField, (field.field_number // 2) % 2) - ] target_phase_even = target_phase target_phase_odd = target_phase @@ -1769,18 +1895,20 @@ def process_chroma( # (field.isFirstField, line_6_burst_present, (field.field_number // 4) % 2) # ] - # offset heterodyne for each line to correct color phase + # this uses the burst measurements to interpolate the correct phase of the color under heterodyne + # phase issues are corrected continiously for each sample using a linear spline interpolated from the burst measurements + # the mixing is performed on the upsampled signal to avoid aliasing introduced from the up-heterodyne mixing product upconvert_chroma_phase_comp( - chroma, - uphet, + chroma, # modifies this in place lineoffset, outwidth, field.phase_sequence, - field.rf.chroma_afc.color_under, - field.rf.chroma_afc.fsc_mhz * 1e6, + field.rf.DecoderParams["color_under_carrier"], + field.rf.SysParams["fsc_mhz"] * 1e6, target_phase_even, - target_phase_odd + target_phase_odd, ) + uphet = chroma else: if field.rf.chroma_afc.conversion_lo is not None: # Explicit conversion LO (ME-SECAM): trim it by the smoothed @@ -1808,6 +1936,7 @@ def process_chroma( else field.rf.chroma_heterodyne ) + uphet = np.zeros((linesout * outwidth), dtype=np.float32) upconvert_chroma( chroma, uphet, @@ -1821,13 +1950,13 @@ def process_chroma( # Mixing the signals will produce waves at the difference and sum of the # frequencies. We only want the difference wave which is at the correct color # carrier frequency here. - # We do however want to be careful to avoid filtering out too much of the sideband. - uphet = sosfiltfilt_rust(field.rf.Filters["FChromaFinal"], uphet) - - # FFT filter way to use a supergauss filter to more sharply cut out the upper harmonic - # This may be a better approach but slows down things a bit much so not using for now - # orig_len = len(uphet) - # uphet = np_fft.irfft(np_fft.rfft(uphet) * field.rf.Filters["FChromaFinal"], n=orig_len) + uphet = filter_chroma_fft( + uphet, + field.rf.SysParams["fsc_mhz"] * 1e6, + field.rf.DecoderParams["color_under_carrier"], + 1.3e6, # lower chroma bandwidth (roughly this for PAL / NTSC) + 80.0 # heterodyne up-mixing attenuation + ) if do_chroma_deemphasis: b, a = field.rf.Filters["chroma_deemphasis"] @@ -1840,19 +1969,121 @@ def process_chroma( else: uphet = comb_c_pal(uphet, outwidth) - # Final automatic chroma gain. - mean_rms = chroma_automatic_gain( + # Chroma AGC + mean_rms, chroma_noise_floor = chroma_automatic_gain( uphet, field.rf.SysParams["burst_abs_ref"], field.phase_sequence, - field.burst_detected_line + field.burst_detected_line, + math.floor(field.usectooutpx(field.rf.SysParams["hsyncPulseUS"])) ) field.rf.field_averages.chroma_level.push(mean_rms) + if field.rf.options.cti_mix != 0: + chroma_transient_improvement( + uphet, + lineoffset * outwidth, + outwidth, + chroma_noise_floor, + field.rf.options.cti_width, + field.rf.options.cti_mix, + ) + return uphet +@njit(cache=True, fastmath=True, nogil=True) +def chroma_transient_improvement( + chroma_data: np.ndarray, + line_start: int, + line_length: int, + base_noise_floor: float, + cti_width: int, + cti_mix: float, +) -> np.ndarray: + """ + Accelerates the sweep rate between color states without warping phase. + Operates symmetrically by measuring forward/backward vector neighbors simultaneously. + """ + # Configure geometric multi-pass decay + decay = 0.25 + num_passes = 4 + + # Pre-calculate mix factors + mix_factors = np.empty(num_passes, dtype=np.float32) + for p in range(num_passes): + mix_factors[p] = cti_mix * (decay ** p) + + # Establish spatial boundaries + remaining_samples = chroma_data.shape[0] - line_start + line_count = remaining_samples // line_length + + # Establish the 4fsc phase-locked sweep radius and scaling threshold + sweep_radius = int(max(4, cti_width * 4)) + mad_threshold = base_noise_floor * math.sqrt(cti_width) + + # Protect against edge bleeding + start_s = sweep_radius + 1 + end_s = line_length - (sweep_radius + 1) + + line_buffer = np.empty(line_length, dtype=chroma_data.dtype) + + # --- Optimized Vector-Ready Loop --- + for l in range(line_count): + curr_line_offset = line_start + (l * line_length) + + for p in range(num_passes): + current_mix = mix_factors[p] + line_buffer[:] = chroma_data[curr_line_offset : curr_line_offset + line_length] + + for s in range(start_s, end_s): + idx = curr_line_offset + s + + i_curr = line_buffer[s] + q_curr = line_buffer[s - 1] + + i_past = line_buffer[s - sweep_radius] + q_past = line_buffer[s - sweep_radius - 1] + + i_future = line_buffer[s + sweep_radius] + q_future = line_buffer[s + sweep_radius - 1] + + # Measure vector distances + i_delta_back = i_curr - i_past + q_delta_back = q_curr - q_past + dist_back = math.sqrt(i_delta_back * i_delta_back + q_delta_back * q_delta_back) + + i_delta_forw = i_future - i_curr + q_delta_forw = q_future - q_curr + dist_forw = math.sqrt(i_delta_forw * i_delta_forw + q_delta_forw * q_delta_forw) + + total_sweep_distance = dist_back + dist_forw + + # Noise Gate + gate_mask = 1.0 if total_sweep_distance > mad_threshold else 0.0 + + # Inherent Boundary (Naturally bounds to [0.0, 1.0) because distances are >= 0) + norm_progress = dist_back / total_sweep_distance if total_sweep_distance != 0 else 0 + + # Transform the progress via sigmoidal sweep-acceleration function + is_lower = norm_progress < 0.5 + + inv_prog = 1.0 - norm_progress + t_low = 4.0 * (norm_progress * norm_progress) + t_high = 1.0 - 4.0 * (inv_prog * inv_prog) + + # Select interpolation weights and anchor to the closer sample (before or after) + t = t_low if is_lower else t_high + anchor_a = i_past if is_lower else i_curr + anchor_b = i_curr if is_lower else i_future + + # 4. Final vector mix + # If gate_mask is 0.0, the delta cancels out and i_curr is cleanly written back + i_target = anchor_a + t * (anchor_b - anchor_a) + chroma_data[idx] = i_curr + (current_mix * gate_mask) * (i_target - i_curr) + + def decode_chroma(field, do_chroma_deemphasis=False): if field.rf.options.write_chroma: """Do track detection if needed and upconvert the chroma signal""" diff --git a/vhsdecode/field.py b/vhsdecode/field.py index e59f0b910..9214d2a7e 100644 --- a/vhsdecode/field.py +++ b/vhsdecode/field.py @@ -1,20 +1,25 @@ import numpy as np +import numba as nb import lddecode.core as ldd import lddecode.utils as lddu from lddecode.utils import inrange from lddecode.utils import hz_to_output_array +import matplotlib.pyplot as plt import vhsdecode.sync as sync import vhsdecode.formats as formats from vhsdecode.doc import detect_dropouts_rf -from vhsdecode.addons.resync import Pulse from vhsdecode.chroma import decode_chroma, decode_chroma_phase_rotation from vhsdecode.debug_plot import plot_data_and_pulses +from collections import namedtuple + NO_PULSES_FOUND = 1 +Pulse = namedtuple('Pulse', ['start', 'len', 'transition', 'level_low', 'level_high']) + # def ynr(data, hpfdata, line_len): # """Dumb vcr-line ynr # """ @@ -1202,22 +1207,6 @@ def refinepulses(self): ) valid_pulses.append((HSYNC, curpulse, good)) i += 1 - elif inrange(curpulse.len, lt_hsync[1], lt_hsync[1] * 3): - # If the pulse is longer than expected, we could have ended up detecting the back - # porch as sync. - # try to move a bit lower to see if we hit a hsync. - data = self.data["video"]["demod_05"][ - curpulse.start : curpulse.start + curpulse.len - ] - threshold = self.rf.iretohz(self.rf.hztoire(data[0]) - 10) - pulses = self.rf.resync.findpulses(data, threshold) - if len(pulses): - newpulse = Pulse(curpulse.start + pulses[0].start, pulses[0].len) - self.rawpulses[i] = newpulse - curpulse = newpulse - else: - # spulse = (HSYNC, self.rawpulses[i], False) - i += 1 elif ( i > 2 and inrange(self.rawpulses[i].len, *lt_eq) @@ -1240,8 +1229,276 @@ def refinepulses(self): return valid_pulses # , num_vblanks - def _try_get_pulses(self, check_levels): - self.rawpulses = self.rf.resync.get_pulses(self, check_levels) + def get_pulses(self, do_level_detect=False): + demod = self.data["video"]["demod_05"] + hsync_len = self.usectoinpx(self.rf.SysParams["hsyncPulseUS"]) + front_porch_len = self.usectoinpx(self.rf.SysParams["activeVideoUS"][0] - self.rf.SysParams["hsyncPulseUS"] - 2) + line_len = round(self.usectoinpx(self.rf.SysParams["line_period"])) + + # 1. Filter out high frequencies + # boxcar FIR filter to remove high frequency data (color burst, pilot tone) + approx_transition = self.usectoinpx(0.22) + window_size = max(3, int(approx_transition)) + if window_size % 2 == 0: + window_size += 1 + + kernel = np.ones(window_size, dtype=np.float64) / window_size + filtered_demod = np.convolve(demod, kernel, mode='same') + + if do_level_detect: + # try to detect the levels by measuring the lower 5%, and 25% of data + n = len(filtered_demod) + idx_5, idx_25 = int(n * 0.05), int(n * 0.25) + partitioned = np.partition(filtered_demod, (idx_5, idx_25)) + sync_tip_est, blanking_est = partitioned[idx_5], partitioned[idx_25] + + pulses, sync_tip_level, blanking_level = FieldShared._get_pulses( + filtered_demod, + hsync_len, + front_porch_len, + line_len, + approx_transition, + sync_tip_est, + blanking_est, + ) + else: + pulses, sync_tip_level, blanking_level = FieldShared._get_pulses( + filtered_demod, + hsync_len, + front_porch_len, + line_len, + approx_transition, + self.rf.DecoderParams["ire0"], + self.rf.DecoderParams["ire0"] + self.rf.DecoderParams["hz_ire"] * self.rf.DecoderParams["vsync_ire"], + ) + + # update levels + if "backporch" in self.rf.options.ire0_adjust: + self.rf.DecoderParams["ire0"] = blanking_level + + if "hsync" in self.rf.options.ire0_adjust: + self.rf.DecoderParams["hz_ire"] = (blanking_level - sync_tip_level) / -self.rf.DecoderParams["vsync_ire"] + + self.sync_tip_level = sync_tip_level + self.blanking_level = blanking_level + + return [Pulse(*p) for p in pulses] + + @staticmethod + @nb.njit(cache=True, nogil=True, fastmath=True) + def _get_pulses( + filtered_demod, + hsync_len, + back_porch_len, + line_len, + approx_transition, + sync_tip_est, + blanking_est, + # tolerance of gap length between hsync pulses (ratio of linelen) + # essentially the speed tolerance in one direction (not +-) + sync_spacing_tolerance=0.15, + min_grid_length=8 # Minimum number of connected lines to keep a grid set + ): + n_samples = len(filtered_demod) + empty_pulses = np.zeros((0, 5), dtype=np.float64) + + slicer_level_est = (sync_tip_est + blanking_est) / 2.0 + + # 2. CANDIDATE EXTRACTION (Time-Domain Width Discriminator) + # Extracts falling/rising edges and ignore noise by enforcing a strict H-sync duration window. + hsync_falls = np.zeros(n_samples // 10, dtype=np.int32) + hsync_rises = np.zeros(n_samples // 10, dtype=np.int32) + cand_count = 0 + w_min, w_max = hsync_len * 0.6, hsync_len * 1.4 + + f_idx = -1 + for i in range(n_samples - 1): + if filtered_demod[i] >= slicer_level_est and filtered_demod[i+1] < slicer_level_est: + f_idx = i + elif f_idx != -1 and filtered_demod[i] < slicer_level_est and filtered_demod[i+1] >= slicer_level_est: + width = i - f_idx + if w_min < width < w_max: + hsync_falls[cand_count] = f_idx + hsync_rises[cand_count] = i + cand_count += 1 + f_idx = -1 + + if cand_count == 0: + return empty_pulses, float(sync_tip_est), float(blanking_est) + + # 3. LOCAL LEVEL MEASUREMENT & AMPLITUDE OUTLIER REJECTION (Clamping Reference Validation) + # Samples local windows to find true sync minima and back porch black levels. + # Uses MAD to purge anomalous pulses caused by VBI equalization lines or static. + cand_sync_levels = np.zeros(cand_count, dtype=np.float64) + cand_porch_levels = np.zeros(cand_count, dtype=np.float64) + + for idx in range(cand_count): + # Extract median sync tip floor to reject impulsive noise + mid = (hsync_falls[idx] + hsync_rises[idx]) // 2 + s_win = filtered_demod[max(0, mid-2) : min(n_samples, mid+3)].copy() + s_win.sort() + cand_sync_levels[idx] = s_win[len(s_win) // 2] if len(s_win) > 0 else sync_tip_est + + # Extract median back porch level for accurate downstream black-level clamping + p_ctr = int(hsync_rises[idx] + back_porch_len * 0.5) + p_win = filtered_demod[max(0, p_ctr-2) : min(n_samples, p_ctr+3)].copy() + p_win.sort() + cand_porch_levels[idx] = p_win[len(p_win) // 2] if len(p_win) > 0 else blanking_est + + # Statistical outlier pruning via Median Absolute Deviation + sync_sort = cand_sync_levels[:cand_count].copy() + sync_sort.sort() + median_sync = sync_sort[cand_count // 2] + + mad_arr = np.abs(cand_sync_levels[:cand_count] - median_sync) + mad_arr.sort() + mad_sync = mad_arr[cand_count // 2] if mad_arr[cand_count // 2] > 0.0 else 1.0 + + amp_count = 0 + for i in range(cand_count): + if abs(cand_sync_levels[i] - median_sync) <= 2.5 * mad_sync: + hsync_falls[amp_count] = hsync_falls[i] + hsync_rises[amp_count] = hsync_rises[i] + cand_sync_levels[amp_count] = cand_sync_levels[i] + cand_porch_levels[amp_count] = cand_porch_levels[i] + amp_count += 1 + + if amp_count == 0: + return empty_pulses, float(sync_tip_est), float(blanking_est) + + # 4. MULTI-GRID DENSITY LOCK (Directional Coherence Filter) + grid_support_count = np.zeros(amp_count, dtype=np.int32) + + # Calculate the directional baseline stride + eff_line_len = line_len * (1.0 + sync_spacing_tolerance) + jitter_tol = line_len * 0.1 + + for i in range(amp_count): + connections = 1 + for j in range(amp_count): + if i == j: + continue + + # Enforce chronological directionality + if j > i: + delta = hsync_falls[j] - hsync_falls[i] + rem = delta % eff_line_len + else: + delta = hsync_falls[i] - hsync_falls[j] + rem = delta % eff_line_len + + # Check if the step falls within the tight jitter window along the directional vector + if (rem < jitter_tol) or (rem > eff_line_len - jitter_tol): + connections += 1 + + grid_support_count[i] = connections + + # Generate a selection mask keeping only elements part of an appropriately sized run + final_mask = np.zeros(amp_count, dtype=np.uint8) + hsync_fit_count = 0 + for i in range(amp_count): + if grid_support_count[i] >= min_grid_length: + final_mask[i] = 1 + hsync_fit_count += 1 + + # 5. REFINED THRESHOLD CALIBRATION + if hsync_fit_count > 0: + sync_sum = 0.0 + porch_sum = 0.0 + for i in range(amp_count): + if final_mask[i]: + sync_sum += cand_sync_levels[i] + porch_sum += cand_porch_levels[i] + + sync_tip_level = float(sync_sum / hsync_fit_count) + back_porch_level = float(porch_sum / hsync_fit_count) + else: + sync_tip_level, back_porch_level = float(sync_tip_est), float(blanking_est) + + # 6. SUBPIXEL SYNTHESIS & COHERENT EDGE EXTRACTION + precise_midpoint = (sync_tip_level + back_porch_level) / 2.0 + f_edges = np.zeros(n_samples // 10, dtype=np.int32) + r_edges = np.zeros(n_samples // 10, dtype=np.int32) + edge_count = 0 + + f_idx = -1 + for i in range(n_samples - 1): + if filtered_demod[i] >= precise_midpoint and filtered_demod[i+1] < precise_midpoint: + f_idx = i + elif f_idx != -1 and filtered_demod[i] < precise_midpoint and filtered_demod[i+1] >= precise_midpoint: + # ------------------------------------------------------------- + # MULTI-GRID DIRECTIONAL NEIGHBORHOOD VALIDATION + # ------------------------------------------------------------- + belongs_to_valid_grid = False + for c_idx in range(amp_count): + if final_mask[c_idx]: + if f_idx >= hsync_falls[c_idx]: + delta = f_idx - hsync_falls[c_idx] + else: + delta = hsync_falls[c_idx] - f_idx + + rem = delta % eff_line_len + if (rem < jitter_tol) or (rem > eff_line_len - jitter_tol): + belongs_to_valid_grid = True + break + + if not belongs_to_valid_grid: + f_idx = -1 + continue + # ------------------------------------------------------------- + f_edges[edge_count], r_edges[edge_count] = f_idx, i + edge_count += 1 + f_idx = -1 + + if edge_count == 0: + return empty_pulses, sync_tip_level, back_porch_level + + # Evaluate the rise-time slope across the 50% slicing point + slopes = np.zeros(edge_count, dtype=np.float64) + sr_count = 0 + for i in range(edge_count): + r = r_edges[i] + if 10 < r < n_samples - 10: + slopes[sr_count] = abs(filtered_demod[r + 1] - filtered_demod[r - 1]) + sr_count += 1 + + if sr_count > 0: + valid_slopes = slopes[:sr_count] + valid_slopes.sort() + fit_sharpness = max(0.1, valid_slopes[sr_count // 2] / max(1e-5, (back_porch_level - sync_tip_level))) + transition = 1.0 / fit_sharpness + else: + transition = approx_transition + + # Map calibrated structures into final subpixel TBC time indices + output_pulses = np.zeros((edge_count, 5), dtype=np.float64) + pulse_count = 0 + for i in range(edge_count): + fe, re = f_edges[i], r_edges[i] + yf0, yf1 = filtered_demod[fe] - precise_midpoint, filtered_demod[fe + 1] - precise_midpoint + subpixel_f = fe + (yf0 / (yf0 - yf1) if (yf0 - yf1) != 0.0 else 0.0) + + yr0, yr1 = filtered_demod[re] - precise_midpoint, filtered_demod[re + 1] - precise_midpoint + subpixel_r = re + (abs(yr0) / (yr1 - yr0) if (yr1 - yr0) != 0.0 else 0.0) + + calc_len = subpixel_r - subpixel_f + if calc_len <= 0: + continue + + s_idx = int(subpixel_f + (calc_len * 0.5)) + p_idx = int(subpixel_r + (back_porch_len * 0.5)) + + output_pulses[pulse_count, 0] = round(subpixel_f) + output_pulses[pulse_count, 1] = round(calc_len) + output_pulses[pulse_count, 2] = transition * 2.0 + output_pulses[pulse_count, 3] = filtered_demod[s_idx] if s_idx < n_samples else sync_tip_level + output_pulses[pulse_count, 4] = filtered_demod[p_idx] if 0 <= p_idx < n_samples else back_porch_level + pulse_count += 1 + + return output_pulses[:pulse_count], sync_tip_level, back_porch_level + + def _try_get_pulses(self, do_level_detect): + self.rawpulses = self.get_pulses(do_level_detect) if ( self.rawpulses is None @@ -1350,16 +1607,127 @@ def _try_get_pulses(self, check_levels): @property def compute_linelocs_issues(self): return self._compute_linelocs_issues + + @staticmethod + @nb.njit(cache=True, fastmath=True, nogil=True) + def _refine_levels_from_vsync_numba(vsync, orig_sync, orig_blank): + # --- Tuning Constants --- + MIN_YIELD_FRACTION = 0.15 + MIN_YIELD_SAMPLES = 5 + MAD_MULTIPLIER = 3.5 + MAD_EPSILON = 1e-5 + AMP_MIN_BOUND = 0.5 + AMP_MAX_BOUND = 1.5 + SIG_POWER_MIN = 1e-5 + VAR_MIN = 1e-6 + MIN_ACCEPTABLE_SNR = 9.0 + SNR_THRESHOLD = 20.0 + + def filter_level_mad(samples, indices): + if samples.size <= MIN_YIELD_SAMPLES: + return samples, indices + + med = np.median(samples) + abs_dev = np.abs(samples - med) + mad = np.median(abs_dev) + + if mad < MAD_EPSILON: + mad = MAD_EPSILON + + clean_mask = abs_dev < (MAD_MULTIPLIER * mad) + return samples[clean_mask], indices[clean_mask] + + total_samples = vsync.size + min_yield = max(int(total_samples * MIN_YIELD_FRACTION), MIN_YIELD_SAMPLES) + + # 1. Hard assignment masks + sync_mask = np.abs(vsync - orig_sync) < np.abs(vsync - orig_blank) + blank_mask = ~sync_mask + indices = np.arange(total_samples) + + # 2. First pass data cleanup using MAD + sync_s, sync_idx = filter_level_mad(vsync[sync_mask], indices[sync_mask]) + blank_s, blank_idx = filter_level_mad(vsync[blank_mask], indices[blank_mask]) + + refined_sync = orig_sync + refined_blank = orig_blank + + # 3. Validation Gate & SNR Calculation + if sync_s.size >= min_yield and blank_s.size >= min_yield: + mean_sync = np.mean(sync_s) + mean_blank = np.mean(blank_s) + + measured_amp = mean_blank - mean_sync + expected_amp = orig_blank - orig_sync + + # Verify transient integrity + if (AMP_MIN_BOUND * expected_amp) < measured_amp < (AMP_MAX_BOUND * expected_amp): + sig_pow = measured_amp ** 2 + if sig_pow < SIG_POWER_MIN: + sig_pow = SIG_POWER_MIN + + # Refine Sync Weight + sync_var = np.var(sync_s) + if sync_var < VAR_MIN: + sync_var = VAR_MIN + + sync_snr = sig_pow / sync_var + if sync_snr >= MIN_ACCEPTABLE_SNR: + sync_weight = sync_snr / (SNR_THRESHOLD + sync_snr) + refined_sync = ((1.0 - sync_weight) * orig_sync) + (sync_weight * mean_sync) + + # Refine Blanking Weight + blank_var = np.var(blank_s) + if blank_var < VAR_MIN: + blank_var = VAR_MIN + + blank_snr = sig_pow / blank_var + if blank_snr >= MIN_ACCEPTABLE_SNR: + blank_weight = blank_snr / (SNR_THRESHOLD + blank_snr) + refined_blank = ((1.0 - blank_weight) * orig_blank) + (blank_weight * mean_blank) + + return refined_sync, refined_blank, sync_s, sync_idx, blank_s, blank_idx + + + def _refine_levels_from_vsync(self, line0loc, meanlinelen): + start = int(round(line0loc + meanlinelen)) + end = int(round(start + meanlinelen * 8.5)) + vsync = self.data["video"]["demod"][start:end] + + orig_sync = self.sync_tip_level + orig_blank = self.blanking_level - def compute_linelocs(self): - has_levels = self.rf.resync.has_levels() + ( + self.sync_tip_level, + self.blanking_level, + sync_s, sync_idx, + blank_s, blank_idx + ) = FieldShared._refine_levels_from_vsync_numba(vsync, orig_sync, orig_blank) - # Skip vsync serration/level detect if we already have levels from a previous field and - # the option is enabled. - # Also run level detection if we encountered issues in this function in the previous field. + if self.rf.debug_plot and self.rf.debug_plot.is_plot_requested("vsync_levels"): + plt.figure(figsize=(11, 5)) + plt.plot(vsync, color='gray', alpha=0.3, label='Raw Signal') + + plt.scatter(sync_idx, sync_s, color='blue', s=2, alpha=0.5, label='Assigned Sync') + plt.scatter(blank_idx, blank_s, color='red', s=2, alpha=0.5, label='Assigned Blanking') + + plt.axhline(y=orig_sync, color='blue', linestyle='--', alpha=0.5, label=f'H Sync ({orig_sync:.3f})') + plt.axhline(y=self.sync_tip_level, color='cyan', linestyle='-', linewidth=2, label=f'V Sync ({self.sync_tip_level:.3f})') + plt.axhline(y=orig_blank, color='red', linestyle='--', alpha=0.5, label=f'H Blanking ({orig_blank:.3f})') + plt.axhline(y=self.blanking_level, color='magenta', linestyle='-', linewidth=2, label=f'V Blanking ({self.blanking_level:.3f})') + + plt.title("Vertical Sync Level Adjustment") + plt.xlabel("Sample") + plt.ylabel("Amplitude") + plt.legend(loc='upper right') + plt.grid(True, alpha=0.2) + plt.tight_layout() + plt.show() + + + def compute_linelocs(self): do_level_detect = ( - not self.rf.options.saved_levels - or not has_levels + self.rf.options.saved_levels is False or self.rf.compute_linelocs_issues is True ) res = self._try_get_pulses(do_level_detect) @@ -1395,7 +1763,6 @@ def compute_linelocs(self): raw_pulses=self.rawpulses, threshold=self.rf.iretohz(self.rf.SysParams["vsync_ire"] / 2), ) - # threshold=self.rf.resync.last_pulse_threshold if first_hsync_loc is None: if self.initphase is False: @@ -1440,6 +1807,21 @@ def compute_linelocs(self): # ldd.logger.info("line0loc %s %s", int(line0loc), int(self.meanlinelen)) + if self.vblank_next is None: + nextfield = linelocs[self.outlinecount - 7] + else: + nextfield = self.vblank_next - (self.inlinelen * 8) + + if not self.rf.options.saved_levels: + # TODO: make into a user facing option to enable / disable vsync levels + self._refine_levels_from_vsync(line0loc, meanlinelen) + + # Apply outputs to decoder parameters + if "backporch" in self.rf.options.ire0_adjust: + self.rf.DecoderParams["ire0"] = self.sync_tip_level + if "hsync" in self.rf.options.ire0_adjust: + self.rf.DecoderParams["hz_ire"] = (self.blanking_level - self.sync_tip_level) / -self.rf.DecoderParams["vsync_ire"] + if self.rf.debug_plot and self.rf.debug_plot.is_plot_requested("line_locs"): line0loc_plot = -1 if line0loc == None else line0loc first_hsync_loc_plot = -1 if first_hsync_loc == None else first_hsync_loc @@ -1463,11 +1845,6 @@ def compute_linelocs(self): extra_lines=[line0loc_plot, first_hsync_loc_plot, vblank_next_plot], ) - if self.vblank_next is None: - nextfield = linelocs[self.outlinecount - 7] - else: - nextfield = self.vblank_next - (self.inlinelen * 8) - if np.count_nonzero(lineloc_errs) < 30: self.rf.compute_linelocs_issues = False elif self.rf.options.saved_levels: @@ -1477,11 +1854,7 @@ def compute_linelocs(self): def refine_linelocs_hsync(self): if not self.rf.options.skip_hsync_refine: - threshold = ( - self.rf.resync.last_pulse_threshold - if self.rf.options.hsync_refine_use_threshold - else self.rf.iretohz(self.rf.SysParams["vsync_ire"] / 2) - ) + threshold = self.rf.iretohz(self.rf.SysParams["vsync_ire"] / 2) return sync.refine_linelocs_hsync(self, self.linebad, threshold) else: @@ -1674,7 +2047,7 @@ def calc_burstmedian(self): def getpulses(self): """Find sync pulses in the demodulated video signal""" - return self.rf.resync.get_pulses(self) + return self.get_pulses() def compute_deriv_error(self, linelocs, baserr): """Disabled this for now as tapes have large variations in line pos diff --git a/vhsdecode/hifi/HiFiDecode.py b/vhsdecode/hifi/HiFiDecode.py index e616dca8b..95d27c759 100644 --- a/vhsdecode/hifi/HiFiDecode.py +++ b/vhsdecode/hifi/HiFiDecode.py @@ -69,20 +69,6 @@ BLOCKS_PER_SECOND = 2 -# assembles the current filter design on a pipe-able filter -class FiltersClass: - def __init__(self, iir_b, iir_a, dtype=REAL_DTYPE): - self.iir_b, self.iir_a = iir_b.astype(dtype), iir_a.astype(dtype) - self.z = lfilter_zi(self.iir_b, self.iir_a) - - def filtfilt(self, data): - return filtfilt(self.iir_b, self.iir_a, data) - - def lfilt(self, data): - output, self.z = lfilter(self.iir_b, self.iir_a, data, zi=self.z) - return output - - @dataclass class AFEParamsVHS: def __init__(self): diff --git a/vhsdecode/main.py b/vhsdecode/main.py index fea9605f8..cbd5607cd 100644 --- a/vhsdecode/main.py +++ b/vhsdecode/main.py @@ -319,6 +319,28 @@ def main(args=None, use_gui=False): default=False, help="Detects and corrects color-under heterodyne rotation change around head-switching area. Corrects chroma artifacts around head-switching area for color-under formats. (Experimental feature)", ) + chroma_group.add_argument( + "--cti_mix", + dest="cti_mix", + type=float, + default=1, + help=( + "Sets Chroma Transient Improvement amount (color-under only). This is wet/dry mix of the overall effect. Set to 0 to disable. Default is 1. Set to 0 to disable CTI." + "\n Chroma Transient Improvement helps to re-focus edges the up-converted color under. This increases how quickly color can change." + ), + ) + chroma_group.add_argument( + "--cti_width", + dest="cti_width", + type=int, + default=2, + help=( + "Sets Chroma Transient Improvement width (color-under only). This controls how sharply to focus the chroma in units of subcarrier cycles. Default is 2." + "\n Since color under has lower bandwidth than the source, some of the hue detail is lost and sharper edges can be recovered with CTI" + "\n * Larger width -> more edge sharpness, less hue detail." + "\n * Smaller width -> more hue detail, less edge sharpness." + ), + ) chroma_group.add_argument( "--dpc", "--disable_phase_correction", @@ -376,7 +398,7 @@ def main(args=None, use_gui=False): " some of the chroma processing." ), ) - plot_options = "demodblock, deemphasis, raw_pulses, line_locs" + plot_options = "demodblock, deemphasis, raw_pulses, line_locs, rf_luma, vsync_levels" debug_group.add_argument( "--dp", "--debug_plot", @@ -393,17 +415,6 @@ def main(args=None, use_gui=False): default=False, help="Disable use of right side of hsync for lineloc detection (old behaviour)", ) - debug_group.add_argument( - "--level_detect_divisor", - dest="level_detect_divisor", - metavar="value", - type=int, - default=3, - help=( - "Use only every nth sample for vsync serration code - may improve speed at" - " cost of minor accuracy. Limited to max 10." - ), - ) debug_group.add_argument( "--no_resample", dest="no_resample", @@ -622,9 +633,10 @@ def build_loader(): rf_options["nldeemp"] = args.nldeemp rf_options["subdeemp"] = args.subdeemp rf_options["y_comb"] = args.y_comb + rf_options["cti_mix"] = args.cti_mix + rf_options["cti_width"] = args.cti_width rf_options["cafc"] = args.cafc rf_options["disable_right_hsync"] = args.disable_right_hsync - rf_options["level_detect_divisor"] = args.level_detect_divisor rf_options["fallback_vsync"] = args.fallback_vsync rf_options["relaxed_line0"] = args.relaxed_line0 rf_options["field_order_confidence"] = int(max(0, min(100, args.field_order_confidence))) diff --git a/vhsdecode/process.py b/vhsdecode/process.py index 6d44ffd83..0dbb34c3c 100644 --- a/vhsdecode/process.py +++ b/vhsdecode/process.py @@ -21,7 +21,6 @@ import vhsdecode.formats as vhs_formats from vhsdecode.addons.chromasep import ChromaSepClass -from vhsdecode.addons.resync import Resync from vhsdecode.addons.chromaAFC import ChromaAFC from vhsdecode.demod import replace_spikes, unwrap_hilbert, smooth_spikes @@ -761,6 +760,8 @@ def __init__( "fm_audio_notch", "chroma_audio_notch", "chroma_offset", + "cti_mix", + "cti_width", "ire0_adjust", "gnrc_afe", "relaxed_line0", @@ -804,6 +805,8 @@ def __init__( rf_options.get("fm_audio_notch", 0) or (tape_format == "HI8"), self.DecoderParams.get("chroma_audio_notch_freq", 0) > 0, int(self.DecoderParams.get("chroma_offset", 5) * (self.freq / 40.0)), + rf_options.get("cti_mix", 1), + rf_options.get("cti_width", 2), ire0_adjust, rf_options.get("gnrc_afe", False), rf_options.get("relaxed_line0", False), @@ -962,32 +965,6 @@ def __init__( # making it through. self.blockcut_end = 1024 - level_detect_divisor = rf_options.get("level_detect_divisor", 1) - - if level_detect_divisor < 1 or level_detect_divisor > 10: - ldd.logger.warning( - "Invalid level detect divisor value %s, using default.", - level_detect_divisor, - ) - level_detect_divisor = 1 - elif inputfreq / level_detect_divisor < 4: - ldd.logger.warning( - "Level detect divisor too high (%s) for input frequency (%s) mhz. Limiting to %s", - level_detect_divisor, - inputfreq, - int(inputfreq // 4), - ) - level_detect_divisor = int(inputfreq // 4) - - self.resync = Resync( - self.freq_hz, - self.SysParams, - self._sysparams_const, - self._processing_thread_pool, - divisor=level_detect_divisor, - debug=self.debug, - ) - if self._chroma_trap: self.chromaTrap = ChromaSepClass( self.freq_hz, self.SysParams["fsc_mhz"], ldd.logger diff --git a/vhsdecode/video_eq.py b/vhsdecode/video_eq.py index 841ab4149..9f45d1d23 100644 --- a/vhsdecode/video_eq.py +++ b/vhsdecode/video_eq.py @@ -17,7 +17,7 @@ def __init__(self, decoder_params, sharpness_level, freq_hz): self._video_eq_filter = { 0: FiltersClass(iir_eq_loband[0], iir_eq_loband[1], freq_hz), - # 1: utils.FiltersClass(iir_eq_hiband[0], iir_eq_hiband[1], freq_hz), + # 1: FiltersClass(iir_eq_hiband[0], iir_eq_hiband[1], freq_hz), } self._gain = decoder_params["video_eq"]["loband"]["order_limit"]