class EndStatistics:
    """Statistics about the 5' or 3' end"""

    def __init__(self, adapter: "SingleAdapter"):
        self.max_error_rate: float = adapter.max_error_rate
        self.sequence: str = adapter.sequence
        self.effective_length: int = adapter.effective_length
        self.has_wildcards: bool = adapter.adapter_wildcards
        self.indels: bool = adapter.indels
        self.adapter_type: str = adapter.descriptive_identifier()
        self.allows_partial_matches: bool = adapter.allows_partial_matches
        # self.errors[l][e] == n iff a sequence of length l matching at e errors was removed n times
        self.errors: Dict[int, Dict[int, int]] = defaultdict(returns_defaultdict_int)
        self.adjacent_bases = {"A": 0, "C": 0, "G": 0, "T": 0, "": 0}
        # TODO avoid hard-coding the list of classes
        self._remove_prefix = isinstance(adapter, FrontAdapter)

    def __repr__(self):
        errors = {k: dict(v) for k, v in self.errors.items()}
        return "EndStatistics(max_error_rate={}, errors={}, adjacent_bases={})".format(
            self.max_error_rate,
            errors,
            self.adjacent_bases,
        )

    def __iadd__(self, other: Any):
        if not isinstance(other, self.__class__):
            raise ValueError("Cannot compare")
        if (
            self.max_error_rate != other.max_error_rate
            or self.sequence != other.sequence
            or self.effective_length != other.effective_length
            or self.indels != other.indels
        ):
            raise RuntimeError("Incompatible EndStatistics, cannot be added")
        for base in ("A", "C", "G", "T", ""):
            self.adjacent_bases[base] += other.adjacent_bases[base]
        for length, error_dict in other.errors.items():
            for errors in error_dict:
                self.errors[length][errors] += other.errors[length][errors]
        return self

    @property
    def lengths(self):
        d = {length: sum(errors.values()) for length, errors in self.errors.items()}
        return d

    def random_match_probabilities(self, gc_content: float) -> List[float]:
        """
        Estimate probabilities that this adapter end matches a
        random sequence. Indels are not taken into account.

        Returns a list p, where p[i] is the probability that
        i bases of this adapter match a random sequence with
        GC content gc_content.
        """
        assert 0.0 <= gc_content <= 1.0
        seq = self.sequence
        # FIXME this is broken for 'anywhere' adapters
        if self._remove_prefix:
            seq = seq[::-1]
        allowed_bases = "CGRYSKMBDHVN" if self.has_wildcards else "GC"
        p = 1.0
        probabilities = [p]
        for i, c in enumerate(seq):
            if c in allowed_bases:
                p *= gc_content / 2.0
            else:
                p *= (1.0 - gc_content) / 2.0
            probabilities.append(p)
        return probabilities class AdapterStatistics(ABC):

    reverse_complemented: int = 0
    name: str
    adapter: "Adapter"

    @abstractmethod
    def __iadd__(self, other):
        pass

    @abstractmethod
    def end_statistics(self) -> Tuple[Optional[EndStatistics], Optional[EndStatistics]]:
        pass

    @abstractmethod
    def add_match(self, match) -> None:
        pass class FrontAdapter(SingleAdapter):
    """A 5' adapter"""

    description = "regular 5'"

    def __init__(self, *args, **kwargs):
        self._force_anywhere = kwargs.pop("force_anywhere", False)
        super().__init__(*args, **kwargs)

    def descriptive_identifier(self) -> str:
        return "regular_five_prime"

    def _aligner(self) -> Aligner:
        return self._make_aligner(
            Where.ANYWHERE.value if self._force_anywhere else Where.FRONT.value
        )

    def match_to(self, sequence: str):
        """
        Attempt to match this adapter to the given read.

        Return a Match instance if a match was found;
        return None if no match was found given the matching criteria (minimum
        overlap length, maximum error rate).
        """
        alignment: Optional[Tuple[int, int, int, int, int, int]] = self.aligner.locate(
            sequence
        )
        if self._debug:
            print(self.aligner.dpmatrix)
        if alignment is None:
            return None
        return RemoveBeforeMatch(*alignment, adapter=self, sequence=sequence)

    def spec(self) -> str:
        return f"{self.sequence}..."

    def create_statistics(self) -> FrontAdapterStatistics:
        return FrontAdapterStatistics(self) class BackAdapter(SingleAdapter):
    """A 3' adapter"""

    description = "regular 3'"

    def __init__(self, *args, **kwargs):
        self._force_anywhere = kwargs.pop("force_anywhere", False)
        super().__init__(*args, **kwargs)

    def descriptive_identifier(self) -> str:
        return "regular_three_prime"

    def _aligner(self):
        return self._make_aligner(
            Where.ANYWHERE.value if self._force_anywhere else Where.BACK.value
        )

    def match_to(self, sequence: str):
        """
        Attempt to match this adapter to the given read.

        Return a Match instance if a match was found;
        return None if no match was found given the matching criteria (minimum
        overlap length, maximum error rate).
        """
        alignment: Optional[Tuple[int, int, int, int, int, int]] = self.aligner.locate(
            sequence
        )
        if self._debug:
            print(self.aligner.dpmatrix)  # pragma: no cover
        if alignment is None:
            return None
        return RemoveAfterMatch(*alignment, adapter=self, sequence=sequence)

    def spec(self) -> str:
        return f"{self.sequence}"

    def create_statistics(self) -> BackAdapterStatistics:
        return BackAdapterStatistics(self) class AnywhereAdapter(SingleAdapter):
    """
    An adapter that can be 5' or 3'. If a match involves the first base of
    the read, it is assumed to be a 5' adapter and a 3' otherwise.
    """

    description = "variable 5'/3'"

    def descriptive_identifier(self) -> str:
        return "anywhere"

    def _aligner(self):
        return self._make_aligner(Where.ANYWHERE.value)

    def match_to(self, sequence: str):
        """
        Attempt to match this adapter to the given string.

        Return a Match instance if a match was found;
        return None if no match was found given the matching criteria (minimum
        overlap length, maximum error rate).
        """
        alignment = self.aligner.locate(sequence.upper())
        if self._debug:
            print(self.aligner.dpmatrix)
        if alignment is None:
            return None
        # guess: if alignment starts at pos 0, it’s a 5' adapter
        if alignment[2] == 0:  # index 2 is rstart
            match = RemoveBeforeMatch(*alignment, adapter=self, sequence=sequence)  # type: ignore
        else:
            match = RemoveAfterMatch(*alignment, adapter=self, sequence=sequence)  # type: ignore
        return match

    def spec(self) -> str:
        return f"...{self.sequence}..."

    def create_statistics(self) -> AnywhereAdapterStatistics:
        return AnywhereAdapterStatistics(self) class LinkedAdapter(Adapter):
    """A 5' adapter combined with a 3' adapter"""

    description = "linked"

    def __init__(
        self,
        front_adapter: SingleAdapter,
        back_adapter: SingleAdapter,
        front_required: bool,
        back_required: bool,
        name: str,
    ):
        super().__init__(name)
        self.front_required = front_required
        self.back_required = back_required

        # The following attributes are needed for the report
        self.where = "linked"
        self.name = _generate_adapter_name() if name is None else name
        self.front_adapter = front_adapter
        self.front_adapter.name = self.name
        self.back_adapter = back_adapter

    def descriptive_identifier(self) -> str:
        return "linked"

    def enable_debug(self):
        self.front_adapter.enable_debug()
        self.back_adapter.enable_debug()

    def match_to(self, sequence: str) -> Optional[LinkedMatch]:
        """
        Match the two linked adapters against a string
        """
        front_match = self.front_adapter.match_to(sequence)
        if self.front_required and front_match is None:
            return None
        if front_match is not None:
            sequence = sequence[front_match.trim_slice()]
        back_match = self.back_adapter.match_to(sequence)
        if back_match is None and (self.back_required or front_match is None):
            return None
        return LinkedMatch(front_match, back_match, self)

    def create_statistics(self) -> LinkedAdapterStatistics:
        return LinkedAdapterStatistics(
            self, front=self.front_adapter, back=self.back_adapter
        )

    @property
    def sequence(self):
        return self.front_adapter.sequence + "..." + self.back_adapter.sequence

    @property
    def remove(self):
        return None

    def spec(self) -> str:
        return f"{self.front_adapter.spec()}...{self.back_adapter.spec()}" class SingleAdapter(Adapter, ABC):
    """
    This class can find a single adapter characterized by sequence, error rate,
    type etc. within reads.

    where --  A Where enum value. This influences where the adapter is allowed to appear within the
        read.

    sequence -- The adapter sequence as string. Will be converted to uppercase.
        Also, Us will be converted to Ts.

    max_errors -- Maximum allowed errors (non-negative float). If the values is less than 1, this is
        interpreted as a rate directly and passed to the aligner. If it is 1 or greater, the value
        is converted to a rate by dividing it by the number of non-N characters in the sequence.

        The error rate is the number of errors in the alignment divided by the length
        of the part of the alignment that matches the adapter.

    minimum_overlap -- Minimum length of the part of the alignment
        that matches the adapter.

    read_wildcards -- Whether IUPAC wildcards in the read are allowed.

    adapter_wildcards -- Whether IUPAC wildcards in the adapter are
        allowed.

    name -- optional name of the adapter. If not provided, the name is set to a
        unique number.
    """

    allows_partial_matches: bool = True

    def __init__(
        self,
        sequence: str,
        max_errors: float = 0.1,
        min_overlap: int = 3,
        read_wildcards: bool = False,
        adapter_wildcards: bool = True,
        name: Optional[str] = None,
        indels: bool = True,
    ):
        self.name: str = _generate_adapter_name() if name is None else name
        super().__init__(self.name)
        self._debug: bool = False
        self.sequence: str = sequence.upper().replace("U", "T")
        if not self.sequence:
            raise ValueError("Adapter sequence is empty")
        if max_errors >= 1 and self.sequence.count("N") != len(self.sequence):
            max_errors /= len(self.sequence) - self.sequence.count("N")
        self.max_error_rate: float = max_errors
        self.min_overlap: int = min(min_overlap, len(self.sequence))
        iupac = frozenset("ABCDGHKMNRSTUVWXY")
        if adapter_wildcards and not set(self.sequence) <= iupac:
            for c in self.sequence:
                if c not in iupac:
                    if c == "I":
                        extra = (
                            "For inosine, consider using N instead and please comment "
                            "on <https://github.com/marcelm/cutadapt/issues/546>."
                        )
                    else:
                        extra = "Use only characters 'ABCDGHKMNRSTUVWXY'."
                    raise InvalidCharacter(
                        f"Character '{c}' in adapter sequence '{self.sequence}' is "
                        f"not a valid IUPAC code. {extra}"
                    )
        # Optimization: Use non-wildcard matching if only ACGT is used
        self.adapter_wildcards: bool = adapter_wildcards and not set(
            self.sequence
        ) <= set("ACGT")
        self.read_wildcards: bool = read_wildcards
        self.indels: bool = indels
        self.aligner = self._aligner()

    def _make_aligner(self, flags: int) -> Aligner:
        # TODO
        # Indels are suppressed by setting their cost very high, but a different algorithm
        # should be used instead.
        indel_cost = 1 if self.indels else 100000
        return Aligner(
            self.sequence,
            self.max_error_rate,
            flags=flags,
            wildcard_ref=self.adapter_wildcards,
            wildcard_query=self.read_wildcards,
            indel_cost=indel_cost,
            min_overlap=self.min_overlap,
        )

    def __repr__(self):
        return (
            "<{cls}(name={name!r}, sequence={sequence!r}, "
            "max_error_rate={max_error_rate}, min_overlap={min_overlap}, "
            "read_wildcards={read_wildcards}, "
            "adapter_wildcards={adapter_wildcards}, "
            "indels={indels})>".format(cls=self.__class__.__name__, **vars(self))
        )

    @property
    def effective_length(self) -> int:
        return self.aligner.effective_length

    def enable_debug(self) -> None:
        """
        Print out the dynamic programming matrix after matching a read to an
        adapter.
        """
        self._debug = True
        self.aligner.enable_debug()

    @abstractmethod
    def _aligner(self):
        pass

    @abstractmethod
    def match_to(self, sequence: str):
        """
        Attempt to match this adapter to the given string.

        Return a Match instance if a match was found;
        return None if no match was found given the matching criteria (minimum
        overlap length, maximum error rate).
        """

    def __len__(self) -> int:
        return len(self.sequence) class LinkedAdapterStatistics(AdapterStatistics):
    """
    Statistics about sequences removed by a lined adapter.
    """

    def __init__(
        self,
        adapter: "LinkedAdapter",
        front: "SingleAdapter",
        back: "SingleAdapter",
    ):
        self.name = adapter.name
        self.adapter = adapter
        self.front = EndStatistics(front)
        self.back = EndStatistics(back)
        self.reverse_complemented = 0

    def __repr__(self):
        return f"LinkedAdapterStatistics(name={self.name}, front={self.front}, back={self.back})"

    def __iadd__(self, other: "LinkedAdapterStatistics"):
        if not isinstance(other, self.__class__):
            raise ValueError("Cannot iadd")
        self.front += other.front
        self.back += other.back
        self.reverse_complemented += other.reverse_complemented
        return self

    def add_match(self, match: "LinkedMatch"):
        # TODO this is duplicated code
        if match.front_match:
            self.front.errors[match.front_match.removed_sequence_length()][
                match.errors
            ] += 1
        if match.back_match:
            adjacent_base = match.back_match.adjacent_base()
            self.back.errors[match.back_match.removed_sequence_length()][
                match.errors
            ] += 1
            try:
                self.back.adjacent_bases[adjacent_base] += 1
            except KeyError:
                self.back.adjacent_bases[""] += 1

    def end_statistics(self) -> Tuple[Optional[EndStatistics], Optional[EndStatistics]]:
        return self.front, self.back class FrontAdapterStatistics(SingleAdapterStatistics):
    def add_match(self, match: "RemoveBeforeMatch"):
        self.end.errors[match.removed_sequence_length()][match.errors] += 1

    def end_statistics(self) -> Tuple[Optional[EndStatistics], Optional[EndStatistics]]:
        return self.end, None class BackAdapterStatistics(SingleAdapterStatistics):
    def add_match(self, match: "RemoveAfterMatch"):
        adjacent_base = match.adjacent_base()
        self.end.errors[match.removed_sequence_length()][match.errors] += 1
        try:
            self.end.adjacent_bases[adjacent_base] += 1
        except KeyError:
            self.end.adjacent_bases[""] += 1

    def end_statistics(self) -> Tuple[Optional[EndStatistics], Optional[EndStatistics]]:
        return None, self.end class AnywhereAdapterStatistics(AdapterStatistics):
    """
    Statistics about sequences removed by a lined adapter.
    """

    def __init__(self, adapter: "AnywhereAdapter"):
        self.name = adapter.name
        self.adapter = adapter
        self.front = EndStatistics(adapter)
        self.back = EndStatistics(adapter)
        self.reverse_complemented = 0

    def __repr__(self):
        return f"AnywhereAdapterStatistics(name={self.name}, front={self.front}, back={self.back})"

    def __iadd__(self, other: "AnywhereAdapterStatistics"):
        if not isinstance(other, AnywhereAdapterStatistics):
            raise ValueError("Cannot add")
        self.front += other.front
        self.back += other.back
        self.reverse_complemented += other.reverse_complemented
        return self

    def add_match(self, match: Union["RemoveBeforeMatch", "RemoveAfterMatch"]) -> None:
        # TODO contains duplicated code from the other add_match() methods
        if isinstance(match, RemoveBeforeMatch):
            self.front.errors[match.removed_sequence_length()][match.errors] += 1
        else:
            adjacent_base = match.adjacent_base()
            self.back.errors[match.removed_sequence_length()][match.errors] += 1
            try:
                self.back.adjacent_bases[adjacent_base] += 1
            except KeyError:
                self.back.adjacent_bases[""] += 1

    def end_statistics(self) -> Tuple[Optional[EndStatistics], Optional[EndStatistics]]:
        return self.front, self.back class OneLine:
    """Wrap any value in this class to print it on one line in the JSON file"""

    def __init__(self, value):
        self.value = value class QualityTrimmer(SingleEndModifier):
    def __init__(self, cutoff_front, cutoff_back, base):
        self.cutoff_front = cutoff_front
        self.cutoff_back = cutoff_back
        self.base = base
        self.trimmed_bases = 0

    def __repr__(self):
        return (
            f"QualityTrimmer(cutoff_front={self.cutoff_front}, "
            f"cutoff_back={self.cutoff_back}, base={self.base})"
        )

    def __call__(self, read, info: ModificationInfo):
        start, stop = quality_trim_index(
            read.qualities, self.cutoff_front, self.cutoff_back, self.base
        )
        self.trimmed_bases += len(read) - (stop - start)
        return read[start:stop] class NextseqQualityTrimmer(SingleEndModifier):
    def __init__(self, cutoff, base):
        self.cutoff = cutoff
        self.base = base
        self.trimmed_bases = 0

    def __repr__(self):
        return f"NextseqQualityTrimmer(cutoff={self.cutoff}, base={self.base})"

    def __call__(self, read, info: ModificationInfo):
        stop = nextseq_trim_index(read, self.cutoff, self.base)
        self.trimmed_bases += len(read) - stop
        return read[:stop] class AdapterCutter(SingleEndModifier):
    """
    Repeatedly find one of multiple adapters in reads.
    The number of times the search is repeated is specified by the
    times parameter.
    """

    def __init__(
        self,
        adapters: Sequence[Adapter],
        times: int = 1,
        action: Optional[str] = "trim",
        index: bool = True,
    ):
        """
        action -- What to do with a found adapter:
          None: Do nothing, only update the ModificationInfo appropriately
          "trim": Remove the adapter and down- or upstream sequence depending on adapter type
          "mask": Replace the part of the sequence that would have been removed with "N" bases
          "lowercase": Convert the part of the sequence that would have been removed to lowercase
          "retain": Like "trim", but leave the adapter sequence itself in the read

        index -- if True, an adapter index (for multiple adapters) is created if possible
        """
        self.times = times
        assert action in ("trim", "mask", "lowercase", "retain", None)
        self.action = action
        self.with_adapters = 0
        self.adapter_statistics = {a: a.create_statistics() for a in adapters}
        if index:
            self.adapters = MultipleAdapters(
                self._regroup_into_indexed_adapters(adapters)
            )
        else:
            self.adapters = MultipleAdapters(adapters)
        if action == "retain" and times > 1:
            raise ValueError("'retain' cannot be combined with times > 1")
        if self.times == 1 and self.action == "trim":
            self.match_and_trim = self._match_and_trim_once_action_trim  # type: ignore

    def __repr__(self):
        return (
            "AdapterCutter("
            f"adapters={self.adapters!r}, times={self.times}, action='{self.action}')"
        )

    def _regroup_into_indexed_adapters(self, adapters):
        prefix, suffix, single = self._split_adapters(adapters)
        # For somewhat better backwards compatibility, avoid re-ordering
        # the adapters when we don’t need to
        if len(prefix) > INDEXING_THRESHOLD or len(suffix) > INDEXING_THRESHOLD:
            result = single
            if len(prefix) > 1:
                result.append(IndexedPrefixAdapters(prefix))
            else:
                result.extend(prefix)
            if len(suffix) > 1:
                result.append(IndexedSuffixAdapters(suffix))
            else:
                result.extend(suffix)
            return result
        else:
            return adapters

    @staticmethod
    def _split_adapters(
        adapters: Sequence[SingleAdapter],
    ) -> Tuple[
        Sequence[SingleAdapter], Sequence[SingleAdapter], Sequence[SingleAdapter]
    ]:
        """
        Split adapters into three different categories so that they can possibly be used
        with a MultiAdapter. Return a tuple (prefix, suffix, other), where
        - prefix is a list of all anchored 5' adapters that MultiAdapter would accept
        - suffix is a list of all anchored 3' adapters that MultiAdapter would accept
        - other is a list of all remaining adapters.
        """
        prefix = []  # type: List[SingleAdapter]
        suffix = []  # type: List[SingleAdapter]
        other = []  # type: List[SingleAdapter]
        for a in adapters:
            if IndexedPrefixAdapters.is_acceptable(a):
                prefix.append(a)
            elif IndexedSuffixAdapters.is_acceptable(a):
                suffix.append(a)
            else:
                other.append(a)
        return prefix, suffix, other

    @staticmethod
    def trim_but_retain_adapter(read, matches: Sequence[Match]):
        start, stop = matches[-1].retained_adapter_interval()
        return read[start:stop]

    @staticmethod
    def masked_read(read, matches: Sequence[Match]):
        start, stop = remainder(matches)
        result = read[:]
        result.sequence = (
            "N" * start + read.sequence[start:stop] + "N" * (len(read) - stop)
        )
        return result

    @staticmethod
    def lowercased_read(read, matches: Sequence[Match]):
        start, stop = remainder(matches)
        result = read[:]
        result.sequence = (
            read.sequence[:start].lower()
            + read.sequence[start:stop].upper()
            + read.sequence[stop:].lower()
        )
        return result

    def __call__(self, read, info: ModificationInfo):
        trimmed_read, matches = self.match_and_trim(read)
        if matches:
            self.with_adapters += 1
            for match in matches:
                self.adapter_statistics[match.adapter].add_match(match)
        info.matches.extend(matches)  # TODO extend or overwrite?
        return trimmed_read

    def match_and_trim(self, read):
        """
        Search for the best-matching adapter in a read, perform the requested action
        ('trim', 'mask' etc. as determined by self.action) and return the
        (possibly) modified read.

        *self.times* adapter removal rounds are done. During each round,
        only the best-matching adapter is trimmed. If no adapter was found in a round,
        no further rounds are attempted.

        Return a pair (trimmed_read, matches), where matches is a list of Match instances.
        """
        matches = []
        if self.action == "lowercase":  # TODO this should not be needed
            read.sequence = read.sequence.upper()
        trimmed_read = read
        for _ in range(self.times):
            match = self.adapters.match_to(trimmed_read.sequence)
            if match is None:
                # if nothing found, attempt no further rounds
                break
            matches.append(match)
            trimmed_read = match.trimmed(trimmed_read)

        if not matches:
            return trimmed_read, []

        if self.action == "trim":
            # read is already trimmed, nothing to do
            pass
        elif self.action == "retain":
            trimmed_read = self.trim_but_retain_adapter(read, matches)
        elif self.action == "mask":
            trimmed_read = self.masked_read(read, matches)
        elif self.action == "lowercase":
            trimmed_read = self.lowercased_read(read, matches)
            assert len(trimmed_read.sequence) == len(read)
        elif self.action is None:
            trimmed_read = read[:]

        return trimmed_read, matches

    def _match_and_trim_once_action_trim(self, read):
        """
        Specalization of match_and_trim for the case that self.times == 1 and self.action == 'trim'
        """
        match = self.adapters.match_to(read.sequence)
        if match is not None:
            return match.trimmed(read), [match]
        else:
            return read, [] class PairedAdapterCutter(PairedEndModifier):
    """
    Trim adapters in pairs from R1 and R2.
    """

    def __init__(self, adapters1, adapters2, action="trim"):
        """
        adapters1 -- list of Adapters to be removed from R1
        adapters2 -- list of Adapters to be removed from R2

        Both lists must have the same, non-zero length.
         read pair is trimmed if adapters1[i] is found in R1 and adapters2[i] in R2.

        action -- What to do with a found adapter: None, 'trim', 'lowercase' or 'mask'
        """
        super().__init__()
        if len(adapters1) != len(adapters2):
            raise PairedAdapterCutterError(
                "The number of adapters to trim from R1 and R2 must be the same. "
                "Given: {} for R1, {} for R2".format(len(adapters1), len(adapters2))
            )
        if not adapters1:
            raise PairedAdapterCutterError("No adapters given")
        self._adapter_pairs = list(zip(adapters1, adapters2))
        logger.debug("Adapter pairs:")
        for a1, a2 in self._adapter_pairs:
            logger.debug(" • %s=%s -- %s=%s", a1.name, a1.spec(), a2.name, a2.spec())
        self.action = action
        self.with_adapters = 0
        self.adapter_statistics = [None, None]
        self.adapter_statistics[0] = {a: a.create_statistics() for a in adapters1}
        self.adapter_statistics[1] = {a: a.create_statistics() for a in adapters2}

    def __repr__(self):
        return f"PairedAdapterCutter(adapter_pairs={self._adapter_pairs!r})"

    def __call__(self, read1, read2, info1, info2):
        """ """
        best_matches = self._find_best_match_pair(read1.sequence, read2.sequence)
        if best_matches is None:
            return read1, read2
        match1, match2 = best_matches
        self.with_adapters += 1
        result = []
        for i, match, read in zip([0, 1], [match1, match2], [read1, read2]):
            trimmed_read = read
            if self.action == "lowercase":
                trimmed_read.sequence = trimmed_read.sequence.upper()

            trimmed_read = match.trimmed(trimmed_read)
            self.adapter_statistics[i][match.adapter].add_match(match)

            if self.action == "trim":
                # read is already trimmed, nothing to do
                pass
            elif self.action == "mask":
                trimmed_read = AdapterCutter.masked_read(read, [match])
            elif self.action == "lowercase":
                trimmed_read = AdapterCutter.lowercased_read(read, [match])
                assert len(trimmed_read.sequence) == len(read)
            elif self.action == "retain":
                trimmed_read = AdapterCutter.trim_but_retain_adapter(read, [match])
            elif self.action is None:  # --no-trim
                trimmed_read = read[:]
            result.append(trimmed_read)
        info1.matches.append(match1)
        info2.matches.append(match2)
        return result

    def _find_best_match_pair(
        self, sequence1: str, sequence2: str
    ) -> Optional[Tuple[Match, Match]]:
        best = None
        best_matches = None
        best_errors = None
        for adapter1, adapter2 in self._adapter_pairs:
            match1 = adapter1.match_to(sequence1)
            if match1 is None:
                continue
            match2 = adapter2.match_to(sequence2)
            if match2 is None:
                continue
            total_matches = match1.matches + match2.matches
            total_errors = match1.errors + match2.errors
            if (
                best is None
                or total_matches > best_matches
                or (total_matches == best_matches and total_errors < best_errors)
            ):
                best = match1, match2
                best_matches = total_matches
                best_errors = total_errors
        return best class ReverseComplementer(SingleEndModifier):
    """Trim adapters from a read and its reverse complement"""

    def __init__(self, adapter_cutter: AdapterCutter, rc_suffix: Optional[str] = " rc"):
        """
        rc_suffix -- suffix to add to the read name if sequence was reverse-complemented
        """
        self.adapter_cutter = adapter_cutter
        self.reverse_complemented = 0
        self._suffix = rc_suffix

    def __repr__(self):
        return f"ReverseComplementer(adapter_cutter={self.adapter_cutter})"

    def __call__(self, read, info: ModificationInfo):
        reverse_read = reverse_complemented_sequence(read)

        forward_trimmed_read, forward_matches = self.adapter_cutter.match_and_trim(read)
        reverse_trimmed_read, reverse_matches = self.adapter_cutter.match_and_trim(
            reverse_read
        )

        forward_match_count = sum(m.matches for m in forward_matches)
        reverse_match_count = sum(m.matches for m in reverse_matches)
        use_reverse_complement = reverse_match_count > forward_match_count

        if use_reverse_complement:
            self.reverse_complemented += 1
            assert reverse_matches
            trimmed_read, matches = reverse_trimmed_read, reverse_matches
            info.is_rc = True
            if self._suffix:
                trimmed_read.name += self._suffix
        else:
            info.is_rc = False
            trimmed_read, matches = forward_trimmed_read, forward_matches

        if matches:
            self.adapter_cutter.with_adapters += 1
            for match in matches:
                stats = self.adapter_cutter.adapter_statistics[match.adapter]
                stats.add_match(match)
                stats.reverse_complemented += bool(use_reverse_complement)
            info.matches.extend(matches)  # TODO extend or overwrite?
        return trimmed_read class PairedEndModifierWrapper(PairedEndModifier):
    """
    Wrapper for modifiers that work on both reads in a paired-end read
    """

    paired = True

    def __init__(
        self,
        modifier1: Optional[SingleEndModifier],
        modifier2: Optional[SingleEndModifier],
    ):
        """Set one of the modifiers to None to work on R1 or R2 only"""
        self._modifier1 = modifier1
        self._modifier2 = modifier2
        if self._modifier1 is None and self._modifier2 is None:
            raise ValueError("Not both modifiers may be None")

    def __repr__(self):
        return (
            "PairedEndModifierWrapper("
            f"modifier1={self._modifier1!r}, modifier2={self._modifier2!r})"
        )

    def __call__(self, read1, read2, info1: ModificationInfo, info2: ModificationInfo):
        if self._modifier1 is None:
            return read1, self._modifier2(read2, info2)  # type: ignore
        if self._modifier2 is None:
            return self._modifier1(read1, info1), read2
        return self._modifier1(read1, info1), self._modifier2(read2, info2) class ReadLengthStatistics:
    """
    Keep track of the lengths of written reads or read pairs
    """

    def __init__(self) -> None:
        # It would be more natural to use a Counter, but a
        # defaultdict is much faster
        self._written_lengths1: DefaultDict[int, int] = defaultdict(int)
        self._written_lengths2: DefaultDict[int, int] = defaultdict(int)

    def update(self, read) -> None:
        """Add a single-end read to the statistics"""
        self._written_lengths1[len(read)] += 1

    def update2(self, read1, read2) -> None:
        """Add a paired-end read to the statistics"""
        self._written_lengths1[len(read1)] += 1
        self._written_lengths2[len(read2)] += 1

    def written_reads(self) -> int:
        """Return number of written reads or read pairs"""
        return sum(self._written_lengths1.values())

    def written_bp(self) -> Tuple[int, int]:
        return (
            self._compute_total_bp(self._written_lengths1),
            self._compute_total_bp(self._written_lengths2),
        )

    def written_lengths(self) -> Tuple[Counter, Counter]:
        return (Counter(self._written_lengths1), Counter(self._written_lengths2))

    @staticmethod
    def _compute_total_bp(counts: DefaultDict[int, int]) -> int:
        return sum(length * count for length, count in counts.items())

    def __iadd__(self, other):
        written_lengths1, written_lengths2 = other.written_lengths()
        for length, count in written_lengths1.items():
            self._written_lengths1[length] += count
        for length, count in written_lengths2.items():
            self._written_lengths2[length] += count
        return self class SingleEndFilter(SingleEndStep):
    """
    A pipeline step that can filter reads, can redirect filtered ones to a writer, and
    counts how many were filtered.
    """

    def __init__(self, writer, predicate: Predicate):
        super().__init__()
        self.filtered = 0
        self.writer = writer
        self.predicate = predicate

    def __repr__(self):
        return f"SingleEndFilter(writer={self.writer}, predicate={self.predicate})"

    def descriptive_identifier(self) -> str:
        return self.predicate.descriptive_identifier()

    def __call__(self, read, info: ModificationInfo) -> bool:
        if self.predicate.test(read, info):
            self.filtered += 1
            if self.writer is not None:
                self.writer.write(read)
            return DISCARD
        return KEEP class PairedEndFilter(PairedEndStep):
    """
    A pipeline step that can filter paired-end reads, redirect them to a file, and counts
    how many read pairs were filtered.

    Different filtering styles are supported, differing by which of the
    two reads in a pair have to fulfill the filtering criterion.
    """

    def __init__(
        self,
        writer,
        predicate1: Predicate,
        predicate2: Predicate,
        pair_filter_mode="any",
    ):
        """
        pair_filter_mode -- these values are allowed:
            'any': The pair is discarded if any read matches.
            'both': The pair is discarded if both reads match.
            'first': The pair is discarded if the first read matches.
        """
        super().__init__()
        if pair_filter_mode not in ("any", "both", "first"):
            raise ValueError("pair_filter_mode must be 'any', 'both' or 'first'")
        self._pair_filter_mode = pair_filter_mode
        self.filtered = 0
        self.writer = writer
        self.predicate1 = predicate1
        self.predicate2 = predicate2
        self._is_filtered: Any
        if predicate2 is None:
            self._is_filtered = self._is_filtered_first
        elif predicate1 is None:
            self._is_filtered = self._is_filtered_second
        elif pair_filter_mode == "any":
            self._is_filtered = self._is_filtered_any
        elif pair_filter_mode == "both":
            self._is_filtered = self._is_filtered_both
        else:
            self._is_filtered = self._is_filtered_first

    def __repr__(self):
        return (
            f"PairedEndFilter(writer={self.writer}, predicate1={self.predicate1}, "
            f"predicate2={self.predicate2}, pair_filter_mode='{self._pair_filter_mode}')"
        )

    def descriptive_identifier(self) -> str:
        if self.predicate1 is not None:
            return self.predicate1.descriptive_identifier()
        else:
            return self.predicate2.descriptive_identifier()

    def _is_filtered_any(
        self, read1, read2, info1: ModificationInfo, info2: ModificationInfo
    ) -> bool:
        return self.predicate1.test(read1, info1) or self.predicate2.test(read2, info2)

    def _is_filtered_both(
        self, read1, read2, info1: ModificationInfo, info2: ModificationInfo
    ) -> bool:
        return self.predicate1.test(read1, info1) and self.predicate2.test(read2, info2)

    def _is_filtered_first(
        self, read1, read2, info1: ModificationInfo, info2: ModificationInfo
    ) -> bool:
        return self.predicate1.test(read1, info1)

    def _is_filtered_second(
        self, read1, read2, info1: ModificationInfo, info2: ModificationInfo
    ) -> bool:
        return self.predicate2.test(read2, info2)

    def __call__(
        self, read1, read2, info1: ModificationInfo, info2: ModificationInfo
    ) -> bool:
        if self._is_filtered(read1, read2, info1, info2):
            self.filtered += 1
            if self.writer is not None:
                self.writer.write(read1, read2)
            return DISCARD
        return KEEP class HasStatistics(ABC):
    """
    Used for the final steps (sinks), which also need to keep
    track of read length statistics
    """

    @abstractmethod
    def get_statistics(self) -> ReadLengthStatistics:
        pass
import sys
import textwrap
from dataclasses import dataclass
from io import StringIO
from collections import defaultdict
from typing import Any, Optional, List, Dict, Iterator, Tuple
from .adapters import (
    EndStatistics,
    AdapterStatistics,
    FrontAdapter,
    BackAdapter,
    AnywhereAdapter,
    LinkedAdapter,
    SingleAdapter,
    LinkedAdapterStatistics,
    FrontAdapterStatistics,
    BackAdapterStatistics,
    AnywhereAdapterStatistics,
)
from .json import OneLine
from .modifiers import (
    QualityTrimmer,
    NextseqQualityTrimmer,
    AdapterCutter,
    PairedAdapterCutter,
    ReverseComplementer,
    PairedEndModifierWrapper,
)
from .statistics import ReadLengthStatistics
from .steps import SingleEndFilter, PairedEndFilter, HasStatistics
                f"{stats.with_adapters[0]:13,d} ({stats.with_adapters_fraction[0]:.1%})\n"
            )

    if stats.reverse_complemented is not None:
        report += (
            "Reverse-complemented:            "
            "{o.reverse_complemented:13,d} ({o.reverse_complemented_fraction:.1%})\n"
        )

    filter_report = format_filter_report(stats)
    if filter_report:
        report += "\n== Read fate breakdown ==\n"
        report += filter_report

    report += textwrap.dedent(
        """\
    {pairs_or_reads} written (passing filters): {o.written:13,d} ({o.written_fraction:.1%})

    Total basepairs processed: {o.total:13,d} bp
    """
    )
    if stats.paired:
        report += "  Read 1: {o.total_bp[0]:13,d} bp\n"
        report += "  Read 2: {o.total_bp[1]:13,d} bp\n"

    if stats.quality_trimmed is not None:
        report += (
            "Quality-trimmed:           "
            f"{stats.quality_trimmed:13,d} bp ({stats.quality_trimmed_fraction:.1%})\n"
        )
        if stats.paired:
            for i in (0, 1):
                if stats.quality_trimmed_bp[i] is not None:
                    report += f"  Read {i + 1}: {stats.quality_trimmed_bp[i]:13,d} bp\n"

    report += (
        "Total written (filtered):  "
        "{o.total_written_bp:13,d} bp ({o.total_written_bp_fraction:.1%})\n"
    )
    if stats.paired:
        report += "  Read 1: {o.written_bp[0]:13,d} bp\n"
        report += "  Read 2: {o.written_bp[1]:13,d} bp\n"
    pairs_or_reads = "Pairs" if stats.paired else "Reads"
    report = report.format(o=stats, pairs_or_reads=pairs_or_reads)
    print_s(report)

    warning = False
    for which_in_pair in (0, 1):
        for adapter_statistics in stats.adapter_stats[which_in_pair]:
            end_statistics = adapter_statistics.end_statistics()
            if end_statistics[0] is not None:
                total_front = sum(end_statistics[0].lengths.values())
            else:
                total_front = 0
            if end_statistics[1] is not None:
                total_back = sum(end_statistics[1].lengths.values())
            else:
                total_back = 0
            total = total_front + total_back
            reverse_complemented = adapter_statistics.reverse_complemented
            adapter = adapter_statistics.adapter
            if isinstance(adapter, BackAdapter):
                assert total_front == 0
            if isinstance(adapter, FrontAdapter):
                assert total_back == 0

            if stats.paired:
                extra = "First read: " if which_in_pair == 0 else "Second read: "
            else:
                extra = ""

            print_s("=" * 3, extra + "Adapter", adapter_statistics.name, "=" * 3)
            print_s()

            if isinstance(adapter_statistics, LinkedAdapterStatistics):
                print_s(
                    "Sequence: {}...{}; Type: linked; Length: {}+{}; "
                    "5' trimmed: {} times; 3' trimmed: {} times".format(
                        adapter_statistics.front.sequence,
                        adapter_statistics.back.sequence,
                        len(adapter_statistics.front.sequence),
                        len(adapter_statistics.back.sequence),
                        total_front,
                        total_back,
                    ),
                    end="",
                )
            else:
                assert isinstance(adapter, (SingleAdapter, AnywhereAdapter))
                print_s(
                    "Sequence: {}; Type: {}; Length: {}; Trimmed: {} times".format(
                        adapter.sequence,
                        adapter.description,
                        len(adapter.sequence),
                        total,
                    ),
                    end="",
                )
            if stats.reverse_complemented is not None:
                print_s(f"; Reverse-complemented: {reverse_complemented} times")
            else:
                print_s()
            if total == 0:
                print_s()
                continue
            if isinstance(adapter_statistics, AnywhereAdapterStatistics):
                assert isinstance(adapter, AnywhereAdapter)
                print_s(total_front, "times, it overlapped the 5' end of a read")
                print_s(
                    total_back, "times, it overlapped the 3' end or was within the read"
                )
                print_s()
                print_s("Minimum overlap:", adapter.min_overlap)
                print_s(error_ranges(adapter_statistics.front))
                print_s("Overview of removed sequences (5')")
                print_s(histogram(adapter_statistics.front, stats.n, gc_content))
                print_s()
                print_s("Overview of removed sequences (3' or within)")
                print_s(histogram(adapter_statistics.back, stats.n, gc_content))
            elif isinstance(adapter_statistics, LinkedAdapterStatistics):

                assert isinstance(adapter, LinkedAdapter)