 
     matrix_iterator forwardMatrix;
     matrix_iterator backwardMatrix;
-    double *scalingFactors;
+    matrix_iterator scalingFactors;
 
     double sM, sI, bM, bI;
 

            }
        }

//        printf("matrix index %u size %u\n", index, ret);
        matrix_size_output[index] = ret;
    }
};

template <target_system system>
struct stride_hmm_matrix_index
{
    typename vector<system, uint32>::view matrix_index;

    stride_hmm_matrix_index(typename vector<system, uint32>::view matrix_index)
        : matrix_index(matrix_index)
    { }

    CUDA_HOST_DEVICE void operator() (const uint32 matrix_group_index)
    {
        constexpr uint32 stride = baq_stride<system>::stride;

        auto initial_element = matrix_index[matrix_group_index * stride];
        for(uint32 i = 1; i < stride; i++)
        {
            matrix_index[matrix_group_index * stride + i] = initial_element + i;
        }
    }
};

        const uint32 matrix_offset = matrix_index % baq_stride<system>::stride;
        const uint32 matrix_base = matrix_index - matrix_offset;

        const uint32 scaling_offset = scaling_index % baq_stride<system>::stride;
        const uint32 scaling_base = scaling_index - scaling_offset;

        // set up matrix and scaling factor pointers
        forwardMatrix = matrix_iterator(&ctx.baq.forward[matrix_base + matrix_offset]);
        backwardMatrix = matrix_iterator(&ctx.baq.backward[matrix_base + matrix_offset]);
        scalingFactors = matrix_iterator(&ctx.baq.scaling[scaling_base + scaling_offset]);

        // get the windows for the current read
        const auto& hmm_reference_window = ctx.baq.hmm_reference_windows[read_index];
