Skip to content

Commit

Permalink
wip
Browse files Browse the repository at this point in the history
  • Loading branch information
athos committed Jan 11, 2024
1 parent fc38815 commit 9bb9a6c
Showing 1 changed file with 47 additions and 14 deletions.
61 changes: 47 additions & 14 deletions src/cljam/io/cram.clj
Original file line number Diff line number Diff line change
Expand Up @@ -304,6 +304,44 @@
#(assoc % :qname (String. ^bytes (RN)))
(throw (ex-info "Omitted read names are not supported yet." {}))))

(defn record-seq [{::keys [len]} features]
(let [arr (byte-array len (byte (int \N)))]
(loop [[f & more] features, rpos 0, spos 1]
(if f
(let [pos (long (:pos f))]
(dotimes [i (- pos spos)]
(aset arr (+ spos i) (byte (int \N))))
(case (:code f)
:subst (recur more (inc rpos) (inc pos))

(:insertion :softclip)
(let [^bytes bs (:bases f)]
(System/arraycopy bs 0 arr pos (alength bs))
(recur more rpos (inc pos)))

:deletion (recur more (+ rpos (long (:len f))) pos)
:ref-skip (recur more (inc rpos) pos)
:insert-base (do (aset arr pos (byte (:base f)))
(recur more rpos (inc pos)))
nil))
(String. arr)))))

(defn record-qual [{::keys [len flag]} features {:keys [QS]}]
(if (pos? (bit-and (long flag) 0x01))
(let [bb (bb/allocate-lsb-byte-buffer len)]
(dotimes [_ len]
(.put bb (byte (+ 33 (long (QS))))))
(String. (.array bb)))
(let [arr (byte-array len)]
(run! (fn [{:keys [^long pos] :as f}]
(case (:code f)
(:read-base :score) (aset arr (dec pos) (byte (:qual f)))
:scores (let [^bytes scores (:scores f)]
(System/arraycopy ^bytes scores 0 arr (dec pos) (alength scores)))
nil))
features)
(String. arr))))

(defn features->cigar [^long read-len features]
(loop [[f & more :as fs] features
pos 0
Expand Down Expand Up @@ -334,25 +372,20 @@
(some->> (not-empty (persistent! acc')) cigar/simplify)))))

(defn build-read-feature-decoder
[{:keys [FN MQ QS FC FP BA BS IN SC HC PD DL RS BB QQ]}]
[{:keys [FN MQ QS FC FP BA BS IN SC HC PD DL RS BB QQ] :as decoders}]
(fn [{::keys [len] :as record}]
(loop [i (long (FN))
prev-pos 0
fs (transient [])]
(if (zero? i)
(let [features (persistent! fs)
record' (assoc record
:mapq (MQ)
:cigar (->> (or (features->cigar len features)
[[len \M]])
(map (fn [[n op]] (str n op)))
(apply str)))
qual (when (pos? (bit-and (long (::flag record)) 0x01))
(let [bb (bb/allocate-lsb-byte-buffer len)]
(dotimes [_ len]
(.put bb (byte (+ 33 (long (QS))))))
(.array bb)))]
(cond-> record' qual (assoc :qual (String. ^bytes qual))))
(let [features (persistent! fs)]
(assoc record
:qual (record-qual record features decoders)
:mapq (MQ)
:cigar (->> (or (features->cigar len features)
[[len \M]])
(map (fn [[n op]] (str n op)))
(apply str))))
(let [code (char (FC))
pos (+ prev-pos (long (FP)))
data (case code
Expand Down

0 comments on commit 9bb9a6c

Please sign in to comment.