     {
         const CRQ_index idx = batch.crq_index(read_index);
 
-        for(uint32 i = idx.qual_start; i < idx.qual_len; i++)
+        for(uint32 i = idx.qual_start; i < idx.qual_start + idx.qual_len; i++)
         {
             const uint8 baq_i = ctx.baq.qualities[i];
             if (baq_i == uint8(-1))


        hmm_common::setup(hmm_index);

//        printf("hmm_glocal(l_ref=%d qstart=%d, l_query=%d)\n", referenceLength, queryStart, queryLen);
//        printf("ref = { ");
//        for(int c = 0; c < referenceLength; c++)
//        {
//            printf("%c ", from_nvbio::iupac16_to_char(referenceBases[c]));
//        }
//        printf("\n");
//
//        printf("que = { ");
//        for(int c = 0; c < queryLen; c++)
//        {
//            printf("%c ", from_nvbio::iupac16_to_char(queryBases[c]));
//        }
//        printf("\n");
//
//        printf("_iqual = { % 3d % 3d % 3d % 3d % 3d ... % 3d % 3d % 3d % 3d % 3d }\n",
//                inputQualities[0], inputQualities[1], inputQualities[2], inputQualities[3], inputQualities[4],
//                inputQualities[queryLen - 5], inputQualities[queryLen - 4], inputQualities[queryLen - 3], inputQualities[queryLen - 2], inputQualities[queryLen - 1]);
//        printf("c->bw = %d, bw = %d, l_ref = %d, l_query = %d\n", MAX_BAND_WIDTH, bandWidth, referenceLength, queryLen);

        /*** forward ***/
        // f[0]

            {
                int u, v11, v01, v10;
                double e = calcEpsilon(referenceBases[k-1], qyi, inputQualities[queryStart+i-1]);
//                printf("read %d: referenceBases[%d-1] = %c inputQualities[%d+%d-1] = %d qyi = %c -> e = %.4f\n",
//                        read_index,
//                        k,
//                        from_nvbio::iupac16_to_char(referenceBases[k-1]),
//                        queryStart,
