diff --git a/src/AlignmentMatrixFast.h b/src/AlignmentMatrixFast.h
index 619692a3..33639cf1 100644
--- a/src/AlignmentMatrixFast.h
+++ b/src/AlignmentMatrixFast.h
@@ -85,15 +85,16 @@ class AlignmentMatrixFast {
 	MatrixElement * getElementEdit(int const x, int const y);
 
 
-	//FAST FUNCTIONS
-
+	//SAFE FUNCTIONS
 	inline MatrixElement * getElementUp(int const x, int const y)
 	{
-		//TODO: Make unnecessary w. loop unrolling
-		if (y < 0 || x < 0) {
+		/*if (y < 0 || x < 0) {
+		       return &empty;
+		}*/
+		       
+		if (x < 0) {
 		       return &empty;
 		}
-
 		if (x < lastCorridor.offset || x >= (lastCorridor.offset + lastCorridor.length)) {
 		       return &empty;
 		}
@@ -104,31 +105,53 @@ class AlignmentMatrixFast {
 
 	inline MatrixElement * getElementCurr(int const x, int const y)
 	{
-		if (y < 0 || x < 0) {
+		/*if (y < 0 || x < 0) {
 		       return &empty;
 		}
 
+		
 		if (x < currentCorridor.offset || x >= (currentCorridor.offset + currentCorridor.length)) {
 		       return &empty;
+		}*/
+		if (x < 0) {
+		       return &empty;
 		}
-
 		return currentLine + (x - currentCorridor.offset);
 	}
 
 	inline MatrixElement * getElementEditCurr(int const x, int const y)
 	{
-		if (y < 0 || x < 0) {
+		/*if (y < 0 || x < 0) {
 		       throw "";
 		}
 
+
 		if (x < currentCorridor.offset || x >= (currentCorridor.offset + currentCorridor.length)) {
 		       throw "";
-		}
+		}*/
+
+		return currentLine + (x - currentCorridor.offset);
+	}
 
+	//END SAFE FUNCTIONS
+
+	//FAST FUNCTIONS
+	inline MatrixElement * getElementUpFast(int const x, int const y)
+	{
+		return lastLine + (x - lastCorridor.offset);
+	}
+
+	inline MatrixElement * getElementCurrFast(int const x, int const y)
+	{
+		return currentLine + (x - currentCorridor.offset);
+	}
+
+	inline MatrixElement * getElementEditCurrFast(int const x, int const y)
+	{
 		return currentLine + (x - currentCorridor.offset);
 	}
 
-	inline char * getDirectionCurr(int const x, int const y) {
+	inline char * getDirectionCurrFast(int const x, int const y) {
 		CorridorLine line = corridorLines[y];
 		return directionMatrix + line.offsetInMatrix + (x - line.offset);
 	}
diff --git a/src/ConvexAlignFast.cpp b/src/ConvexAlignFast.cpp
index 11024d87..9e23407e 100644
--- a/src/ConvexAlignFast.cpp
+++ b/src/ConvexAlignFast.cpp
@@ -449,7 +449,7 @@ int ConvexAlignFast::SingleAlign(int const mode, CorridorLine * corridorLines,
 //
 ////	Timer t1;
 ////	t1.ST();
-		AlignmentMatrixFast::Score score = FastfwdFillMatrix(refSeq, qrySeq, fwdResults,
+		AlignmentMatrixFast::Score score = FastUnrolledfwdFillMatrixMaster(refSeq, qrySeq, fwdResults,
 				mode);
 ////	fprintf(stderr, "fill: %f\n", t1.ET());
 //
@@ -718,7 +718,7 @@ AlignmentMatrixFast::Score ConvexAlignFast::FastfwdFillMatrix(char const * const
 
 			AlignmentMatrixFast::MatrixElement * current = matrix->getElementEditCurr(x,y);
 
-			char & currentDirection = *matrix->getDirectionCurr(x, y);
+			char & currentDirection = *matrix->getDirectionCurrFast(x, y);
 
 			if (del_run > 0 && max_cell == left_cell) {
 				current->score = max_cell;
@@ -779,344 +779,36 @@ AlignmentMatrixFast::Score ConvexAlignFast::FastUnrolledfwdFillMatrix(char const
 		char const * const qrySeq, FwdResults & fwdResult, int readId) {
 
 	AlignmentMatrixFast::Score curr_max = -1.0f;
-	
-	/*
-	 * BEGIN UNROLLED Y=0 ITERATION
-	 */
-	matrix->prepareLine(0);
 
-	int xOffset = matrix->getCorridorOffset(0);
-	char const read_char_cache = qrySeq[0];
-
-	int xMax=std::min(xOffset + matrix->getCorridorLength(0),matrix->getWidth());
-
-	int x = std::max(0, xOffset);
-
-	AlignmentMatrixFast::Score diag_score = matrix->empty.score;
-	AlignmentMatrixFast::MatrixElement const & up = matrix->empty;
-	AlignmentMatrixFast::MatrixElement const & left = matrix->empty;
-
-	bool const eq = read_char_cache == refSeq[x];
-	AlignmentMatrixFast::Score const diag_cell = diag_score
-			+ ((eq) ? mat : mis);
-
-	AlignmentMatrixFast::Score up_cell = 0;
-	AlignmentMatrixFast::Score left_cell = 0;
-
-	int ins_run = 0;
-	int del_run = 0;
-
-	if (up.direction == CIGAR_I) {
-		ins_run = up.indelRun;
-		if (up.score == 0) {
-			up_cell = 0;
-		} else {
-			up_cell = up.score
-					+ std::min(gap_ext_min,
-							gap_ext + ins_run * gap_decay);
-		}
-	} else {
-		up_cell = up.score + gap_open_read;
-	}
-
-	if (left.direction == CIGAR_D) {
-		del_run = left.indelRun;
-		if (left.score == 0) {
-			left_cell = 0;
-		} else {
-			left_cell = left.score
-					+ std::min(gap_ext_min,
-							gap_ext + del_run * gap_decay);
-		}
-	} else {
-		left_cell = left.score + gap_open_ref;
-	}
-
-	//find max
-	AlignmentMatrixFast::Score max_cell = 0;
-	max_cell = std::max(left_cell, max_cell);
-	max_cell = std::max(diag_cell, max_cell);
-	max_cell = std::max(up_cell, max_cell);
-
-	AlignmentMatrixFast::MatrixElement * current = matrix->getElementEditCurr(x,0);
-
-	char & currentDirection = *matrix->getDirectionCurr(x, 0);
-
-	if (del_run > 0 && max_cell == left_cell) {
-		current->score = max_cell;
-		current->direction = CIGAR_D;
-		currentDirection = CIGAR_D;
-		current->indelRun = del_run + 1;
-	} else if (ins_run > 0 && max_cell == up_cell) {
-		current->score = max_cell;
-		current->direction = CIGAR_I;
-		currentDirection = CIGAR_I;
-		current->indelRun = ins_run + 1;
-	} else if (max_cell == diag_cell) {
-		current->score = max_cell;
-		if (eq) {
-			current->direction = CIGAR_EQ;
-			currentDirection = CIGAR_EQ;
-		} else {
-			current->direction = CIGAR_X;
-			currentDirection = CIGAR_X;
-		}
-		current->indelRun = 0;
-	} else if (max_cell == left_cell) {
-		current->score = max_cell;
-		current->direction = CIGAR_D;
-		currentDirection = CIGAR_D;
-		current->indelRun = 1;
-	} else if (max_cell == up_cell) {
-		current->score = max_cell;
-		current->direction = CIGAR_I;
-		currentDirection = CIGAR_I;
-		current->indelRun = 1;
-	} else {
-		current->score = 0;
-		current->direction = CIGAR_STOP;
-		currentDirection = CIGAR_STOP;
-		current->indelRun = 0;
-	}
-
-	if (max_cell > curr_max) {
-		curr_max = max_cell;
-		fwdResult.best_ref_index = x;
-		fwdResult.best_read_index = 0;
-		fwdResult.max_score = curr_max;
-	}
-
-
-	x++;
-
-	for (; x < xMax; ++x) {
-
-		AlignmentMatrixFast::Score diag_score = matrix->empty.score;
-		AlignmentMatrixFast::MatrixElement const & up = matrix->empty;
-		AlignmentMatrixFast::MatrixElement const & left = *matrix->getElementCurr(x - 1, 0);
-
-		bool const eq = read_char_cache == refSeq[x];
-		AlignmentMatrixFast::Score const diag_cell = diag_score
-				+ ((eq) ? mat : mis);
-
-		AlignmentMatrixFast::Score up_cell = 0;
-		AlignmentMatrixFast::Score left_cell = 0;
-
-		int ins_run = 0;
-		int del_run = 0;
-
-		if (up.direction == CIGAR_I) {
-			ins_run = up.indelRun;
-			if (up.score == 0) {
-				up_cell = 0;
-			} else {
-				up_cell = up.score
-						+ std::min(gap_ext_min,
-								gap_ext + ins_run * gap_decay);
-			}
-		} else {
-			up_cell = up.score + gap_open_read;
-		}
-
-		if (left.direction == CIGAR_D) {
-			del_run = left.indelRun;
-			if (left.score == 0) {
-				left_cell = 0;
-			} else {
-				left_cell = left.score
-						+ std::min(gap_ext_min,
-								gap_ext + del_run * gap_decay);
-			}
-		} else {
-			left_cell = left.score + gap_open_ref;
-		}
-
-		//find max
-		AlignmentMatrixFast::Score max_cell = 0;
-		max_cell = std::max(left_cell, max_cell);
-		max_cell = std::max(diag_cell, max_cell);
-		max_cell = std::max(up_cell, max_cell);
-
-		AlignmentMatrixFast::MatrixElement * current = matrix->getElementEditCurr(x,0);
-
-		char & currentDirection = *matrix->getDirectionCurr(x, 0);
-
-		if (del_run > 0 && max_cell == left_cell) {
-			current->score = max_cell;
-			current->direction = CIGAR_D;
-			currentDirection = CIGAR_D;
-			current->indelRun = del_run + 1;
-		} else if (ins_run > 0 && max_cell == up_cell) {
-			current->score = max_cell;
-			current->direction = CIGAR_I;
-			currentDirection = CIGAR_I;
-			current->indelRun = ins_run + 1;
-		} else if (max_cell == diag_cell) {
-			current->score = max_cell;
-			if (eq) {
-				current->direction = CIGAR_EQ;
-				currentDirection = CIGAR_EQ;
-			} else {
-				current->direction = CIGAR_X;
-				currentDirection = CIGAR_X;
-			}
-			current->indelRun = 0;
-		} else if (max_cell == left_cell) {
-			current->score = max_cell;
-			current->direction = CIGAR_D;
-			currentDirection = CIGAR_D;
-			current->indelRun = 1;
-		} else if (max_cell == up_cell) {
-			current->score = max_cell;
-			current->direction = CIGAR_I;
-			currentDirection = CIGAR_I;
-			current->indelRun = 1;
-		} else {
-			current->score = 0;
-			current->direction = CIGAR_STOP;
-			currentDirection = CIGAR_STOP;
-			current->indelRun = 0;
-		}
-
-		if (max_cell > curr_max) {
-			curr_max = max_cell;
-			fwdResult.best_ref_index = x;
-			fwdResult.best_read_index = 0;
-			fwdResult.max_score = curr_max;
-		}
-	}
-
-	/*
-	 * END UNROLLED Y=0 ITERATION
-	 */
-
-	/*
-	 * BEGIN MAIN Y LOOP
-	 */
-	for (int y = 1; y < matrix->getHeight(); ++y) {
+	for (int y = 0; y < matrix->getHeight(); ++y) {
 
 		matrix->prepareLine(y);
-
 		int xOffset = matrix->getCorridorOffset(y);
-
 		char const read_char_cache = qrySeq[y];
-
 		int xMax=std::min(xOffset + matrix->getCorridorLength(y),matrix->getWidth());
 
-		int x = std::max(0, xOffset);
-
-		/*
-		 * BEGIN UNROLLED Y=Y, X=0 ITERATION
-		 */
-		AlignmentMatrixFast::Score diag_score = matrix->empty.score;
-		AlignmentMatrixFast::MatrixElement const & up = *matrix->getElementUp(x,y - 1);
-		AlignmentMatrixFast::MatrixElement const & left = matrix->empty;
-
-		bool const eq = read_char_cache == refSeq[x];
-		AlignmentMatrixFast::Score const diag_cell = diag_score
-				+ ((eq) ? mat : mis);
-
-		AlignmentMatrixFast::Score up_cell = 0;
-		AlignmentMatrixFast::Score left_cell = 0;
-
-		int ins_run = 0;
-		int del_run = 0;
-
-		if (up.direction == CIGAR_I) {
-			ins_run = up.indelRun;
-			if (up.score == 0) {
-				up_cell = 0;
-			} else {
-				up_cell = up.score
-						+ std::min(gap_ext_min,
-								gap_ext + ins_run * gap_decay);
-			}
-		} else {
-			up_cell = up.score + gap_open_read;
-		}
-
-		if (left.direction == CIGAR_D) {
-			del_run = left.indelRun;
-			if (left.score == 0) {
-				left_cell = 0;
-			} else {
-				left_cell = left.score
-						+ std::min(gap_ext_min,
-								gap_ext + del_run * gap_decay);
-			}
-		} else {
-			left_cell = left.score + gap_open_ref;
-		}
-
-		//find max
-		AlignmentMatrixFast::Score max_cell = 0;
-		max_cell = std::max(left_cell, max_cell);
-		max_cell = std::max(diag_cell, max_cell);
-		max_cell = std::max(up_cell, max_cell);
-
-		AlignmentMatrixFast::MatrixElement * current = matrix->getElementEditCurr(x,y);
-
-		char & currentDirection = *matrix->getDirectionCurr(x, y);
-
-		if (del_run > 0 && max_cell == left_cell) {
-			current->score = max_cell;
-			current->direction = CIGAR_D;
-			currentDirection = CIGAR_D;
-			current->indelRun = del_run + 1;
-		} else if (ins_run > 0 && max_cell == up_cell) {
-			current->score = max_cell;
-			current->direction = CIGAR_I;
-			currentDirection = CIGAR_I;
-			current->indelRun = ins_run + 1;
-		} else if (max_cell == diag_cell) {
-			current->score = max_cell;
-			if (eq) {
-				current->direction = CIGAR_EQ;
-				currentDirection = CIGAR_EQ;
+		for (int x = std::max(0, xOffset); x < xMax; ++x) {
+			AlignmentMatrixFast::Score diag_score;
+			AlignmentMatrixFast::MatrixElement const * p_up;
+			AlignmentMatrixFast::MatrixElement const * p_left;
+			AlignmentMatrixFast::MatrixElement * current;
+
+			if(y==0 ||
+			  (x-1 <= (matrix->getCorridorOffset(y-1)) ) ||
+			  (x >= (matrix->getCorridorOffset(y-1)+matrix->getCorridorLength(y-1)) )  ) {
+				diag_score = matrix->getElementUp(x - 1, y - 1)->score;
+				p_up = matrix->getElementUp(x,y - 1);
+				p_left = matrix->getElementCurr(x - 1, y);
+				current = matrix->getElementEditCurr(x,y);
 			} else {
-				current->direction = CIGAR_X;
-				currentDirection = CIGAR_X;
+				diag_score = matrix->getElementUpFast(x - 1, y - 1)->score;
+				p_up = matrix->getElementUpFast(x,y - 1);
+				p_left = matrix->getElementCurrFast(x - 1, y);
+				current = matrix->getElementEditCurrFast(x,y);
 			}
-			current->indelRun = 0;
-		} else if (max_cell == left_cell) {
-			current->score = max_cell;
-			current->direction = CIGAR_D;
-			currentDirection = CIGAR_D;
-			current->indelRun = 1;
-		} else if (max_cell == up_cell) {
-			current->score = max_cell;
-			current->direction = CIGAR_I;
-			currentDirection = CIGAR_I;
-			current->indelRun = 1;
-		} else {
-			current->score = 0;
-			current->direction = CIGAR_STOP;
-			currentDirection = CIGAR_STOP;
-			current->indelRun = 0;
-		}
-
-		if (max_cell > curr_max) {
-			curr_max = max_cell;
-			fwdResult.best_ref_index = x;
-			fwdResult.best_read_index = y;
-			fwdResult.max_score = curr_max;
-		}
-
-
-		x++;
-
-		/*
-		 * END UNROLLED Y=Y, X=0 ITERATION
-		 */
-
-		/*
-		 * BEGIN MAIN Y=Y, X=X ITERATION
-		 */
-		for (; x < xMax; ++x) {
-
-			AlignmentMatrixFast::Score diag_score = matrix->getElementUp(x - 1, y - 1)->score;
-			AlignmentMatrixFast::MatrixElement const & up = *matrix->getElementUp(x,y - 1);
-			AlignmentMatrixFast::MatrixElement const & left = *matrix->getElementCurr(x - 1, y);
+			AlignmentMatrixFast::MatrixElement const & up=*p_up;
+			AlignmentMatrixFast::MatrixElement const & left=*p_left;
+			char & currentDirection = *matrix->getDirectionCurrFast(x, y);
 
 			bool const eq = read_char_cache == refSeq[x];
 			AlignmentMatrixFast::Score const diag_cell = diag_score
@@ -1160,10 +852,6 @@ AlignmentMatrixFast::Score ConvexAlignFast::FastUnrolledfwdFillMatrix(char const
 			max_cell = std::max(diag_cell, max_cell);
 			max_cell = std::max(up_cell, max_cell);
 
-			AlignmentMatrixFast::MatrixElement * current = matrix->getElementEditCurr(x,y);
-
-			char & currentDirection = *matrix->getDirectionCurr(x, y);
-
 			if (del_run > 0 && max_cell == left_cell) {
 				current->score = max_cell;
 				current->direction = CIGAR_D;
@@ -1209,9 +897,6 @@ AlignmentMatrixFast::Score ConvexAlignFast::FastUnrolledfwdFillMatrix(char const
 			}
 
 		}
-		/*
-		 * END MAIN Y=Y, X=X ITERATION
-		 */
 
 	}
 	fwdResult.qend = (matrix->getHeight() - fwdResult.best_read_index) - 1;
diff --git a/src/ConvexAlignFast.h b/src/ConvexAlignFast.h
index c3324228..3ecd5e64 100644
--- a/src/ConvexAlignFast.h
+++ b/src/ConvexAlignFast.h
@@ -23,6 +23,7 @@
 #include "IAlignment.h"
 
 #include <stdint.h>
+#include <algorithm>
 
 #include "Types.h"
 #include "AlignmentMatrixFast.h"
@@ -145,6 +146,462 @@ class ConvexAlignFast: public IAlignment {
 	AlignmentMatrixFast::Score FastUnrolledfwdFillMatrix(char const * const refSeq,
 				char const * const qrySeq, FwdResults & fwdResults, int readId);
 
+	/*AlignmentMatrixFast::Score FastUnrolledfwdFillMatrixMaster(char const * const refSeq,
+				char const * const qrySeq, FwdResults & fwdResults, int readId);
+
+	AlignmentMatrixFast::Score FastUnrolledfwdFillMatrixLine(char const * const refSeq,
+				char const * const qrySeq, FwdResults & fwdResults, int readId,
+				AlignmentMatrixFast::Score diag_score,
+				AlignmentMatrixFast::MatrixElement const & up,
+				AlignmentMatrixFast::MatrixElement const & left,
+				AlignmentMatrixFast::MatrixElement * current,
+				char & currentDirection,
+				AlignmentMatrixFast::Score & curr_max,
+				int x,
+				int y);*/
+
+#define FILLMATRIX(DIAG_SCORE,UP_REF,LEFT_REF,CURRENT,CURRENT_DIR) {\
+			diag_score=DIAG_SCORE;\
+			AlignmentMatrixFast::MatrixElement const & up = UP_REF;\
+			AlignmentMatrixFast::MatrixElement const & left = LEFT_REF;\
+			AlignmentMatrixFast::MatrixElement * current=CURRENT;\
+			char & currentDirection=CURRENT_DIR;\
+			bool const eq = read_char_cache == refSeq[x]; \
+			diag_cell = diag_score + (eq ? mat : mis);\
+\
+			up_cell = 0;\
+			left_cell = 0;\
+\
+			ins_run = 0;\
+			del_run = 0;\
+\
+			if (up.direction == CIGAR_I) {\
+				ins_run = up.indelRun;\
+				if (up.score == 0) {\
+					up_cell = 0;\
+				} else {\
+					up_cell = up.score\
+							+ std::min(gap_ext_min,\
+									gap_ext + ins_run * gap_decay);\
+				}\
+			} else {\
+				up_cell = up.score + gap_open_read;\
+			}\
+\
+			if (left.direction == CIGAR_D) {\
+				del_run = left.indelRun;\
+				if (left.score == 0) {\
+					left_cell = 0;\
+				} else {\
+					left_cell = left.score\
+							+ std::min(gap_ext_min,\
+									gap_ext + del_run * gap_decay);\
+				}\
+			} else {\
+				left_cell = left.score + gap_open_ref;\
+			}\
+\
+			AlignmentMatrixFast::Score max_cell = 0;\
+			max_cell = std::max(left_cell, max_cell);\
+			max_cell = std::max(diag_cell, max_cell);\
+			max_cell = std::max(up_cell, max_cell);\
+\
+			if (del_run > 0 && max_cell == left_cell) {\
+				current->score = max_cell;\
+				current->direction = CIGAR_D;\
+				currentDirection = CIGAR_D;\
+				current->indelRun = del_run + 1;\
+			} else if (ins_run > 0 && max_cell == up_cell) {\
+				current->score = max_cell;\
+				current->direction = CIGAR_I;\
+				currentDirection = CIGAR_I;\
+				current->indelRun = ins_run + 1;\
+			} else if (max_cell == diag_cell) {\
+				current->score = max_cell;\
+				if (eq) {\
+					current->direction = CIGAR_EQ;\
+					currentDirection = CIGAR_EQ;\
+				} else {\
+					current->direction = CIGAR_X;\
+					currentDirection = CIGAR_X;\
+				}\
+				current->indelRun = 0;\
+			} else if (max_cell == left_cell) {\
+				current->score = max_cell;\
+				current->direction = CIGAR_D;\
+				currentDirection = CIGAR_D;\
+				current->indelRun = 1;\
+			} else if (max_cell == up_cell) {\
+				current->score = max_cell;\
+				current->direction = CIGAR_I;\
+				currentDirection = CIGAR_I;\
+				current->indelRun = 1;\
+			} else {\
+				current->score = 0;\
+				current->direction = CIGAR_STOP;\
+				currentDirection = CIGAR_STOP;\
+				current->indelRun = 0;\
+			}\
+\
+			if (max_cell > curr_max) {\
+				curr_max = max_cell;\
+				fwdResult.best_ref_index = x;\
+				fwdResult.best_read_index = y;\
+				fwdResult.max_score = curr_max;\
+			}\
+		}
+
+	inline AlignmentMatrixFast::Score FastUnrolledfwdFillMatrixMaster(char const * const refSeq,
+			char const * const qrySeq, FwdResults & fwdResult, int readId) {
+
+		AlignmentMatrixFast::Score curr_max = -1.0f;
+		int const matrix_height=matrix->getHeight();
+
+
+		AlignmentMatrixFast::Score diag_score;
+		AlignmentMatrixFast::Score diag_cell;
+		AlignmentMatrixFast::Score left_cell;
+		AlignmentMatrixFast::Score up_cell;
+
+		int ins_run;
+		int del_run;
+
+
+		{
+			int y=0;
+
+			matrix->prepareLine(y);
+			int xOffset = matrix->getCorridorOffset(y);
+			int xMax=std::min(xOffset + matrix->getCorridorLength(y),matrix->getWidth());
+
+			char const read_char_cache = qrySeq[y];
+
+			for (int x = std::max(0, xOffset); x < xMax; ++x) {
+				FILLMATRIX(matrix->empty.score,matrix->empty,*matrix->getElementCurr(x - 1, y),matrix->getElementEditCurrFast(x,y),*matrix->getDirectionCurrFast(x, y));
+			}
+		}
+
+
+		for (int y = 1; y < matrix_height; ++y) {
+
+			matrix->prepareLine(y);
+			int xOffset = matrix->getCorridorOffset(y);
+
+			int xMin=std::max(0, xOffset);
+			int xMax=std::min(xOffset + matrix->getCorridorLength(y),matrix->getWidth());
+
+			int xFastMin=std::max(0,matrix->getCorridorOffset(y-1)+1);
+			int xFastMax=std::min(matrix->getCorridorOffset(y-1)+matrix->getCorridorLength(y-1),matrix->getWidth());
+
+
+			char const read_char_cache = qrySeq[y];
+
+			if(xMin>xFastMin || xFastMax>xMax)
+			{
+				for(int x=xMin; x<xMax; x++)
+				{
+					FILLMATRIX(matrix->getElementUp(x - 1, y - 1)->score,*matrix->getElementUp(x,y - 1),*matrix->getElementCurr(x - 1, y),matrix->getElementEditCurrFast(x,y),*matrix->getDirectionCurrFast(x, y));
+				}			
+			} else {
+				for(int x=xMin; x<=xFastMin; x++)
+				{
+					FILLMATRIX(matrix->getElementUp(x - 1, y - 1)->score,*matrix->getElementUp(x,y - 1),*matrix->getElementCurr(x - 1, y),matrix->getElementEditCurrFast(x,y),*matrix->getDirectionCurrFast(x, y));
+				}
+
+				for(int x=xFastMin+1; x<xFastMax; x++)
+				{
+					FILLMATRIX(matrix->getElementUpFast(x - 1, y - 1)->score,*matrix->getElementUpFast(x,y - 1),*matrix->getElementCurrFast(x - 1, y),matrix->getElementEditCurrFast(x,y),*matrix->getDirectionCurrFast(x, y));
+				}
+
+				for(int x=xFastMax; x<xMax; x++)
+				{
+					FILLMATRIX(matrix->getElementUp(x - 1, y - 1)->score,*matrix->getElementUp(x,y - 1),*matrix->getElementCurr(x - 1, y),matrix->getElementEditCurrFast(x,y),*matrix->getDirectionCurrFast(x, y));
+				}			
+			}
+		}
+		fwdResult.qend = (matrix->getHeight() - fwdResult.best_read_index) - 1;
+		if (matrix->getHeight() == 0) {
+			fwdResult.best_read_index = fwdResult.best_ref_index = 0;
+		}
+
+		return curr_max;
+	}
+
+	// inline AlignmentMatrixFast::Score FastUnrolledfwdFillMatrixMaster(char const * const refSeq,
+	// 		char const * const qrySeq, FwdResults & fwdResult, int readId) {
+
+	// 	AlignmentMatrixFast::Score curr_max = -1.0f;
+	// 	int matrix_height=matrix->getHeight();
+
+
+	// 	{
+	// 		int y=0;
+
+	// 		matrix->prepareLine(y);
+	// 		int xOffset = matrix->getCorridorOffset(y);
+	// 		int xMax=std::min(xOffset + matrix->getCorridorLength(y),matrix->getWidth());
+
+	// 		for (int x = std::max(0, xOffset); x < xMax; ++x) {
+	// 			ConvexAlignFast::FastUnrolledfwdFillMatrixLine(refSeq,qrySeq,fwdResult,readId,
+	// 				matrix->empty.score,
+	// 				matrix->empty,
+	// 				*matrix->getElementCurr(x - 1, y),
+	// 				matrix->getElementEditCurrFast(x,y),
+	// 				*matrix->getDirectionCurrFast(x, y),
+	// 				curr_max,
+	// 				x,
+	// 				y,
+	// 				qrySeq[0]);
+	// 		}
+	// 	}
+
+
+	// 	for (int y = 1; y < matrix_height; ++y) {
+
+	// 		matrix->prepareLine(y);
+	// 		int xOffset = matrix->getCorridorOffset(y);
+
+	// 		int xMin=std::max(0, xOffset);
+	// 		int xMax=std::min(xOffset + matrix->getCorridorLength(y),matrix->getWidth());
+
+	// 		int xFastMin=std::max(0,matrix->getCorridorOffset(y-1)+1);
+	// 		int xFastMax=std::min(matrix->getCorridorOffset(y-1)+matrix->getCorridorLength(y-1),matrix->getWidth());
+
+
+	// 		char const read_char_cache = qrySeq[y];
+
+	// 		if(xMin>xFastMin || xFastMax>xMax)
+	// 		{
+	// 			for(int x=xMin; x<xMax; x++)
+	// 			{
+	// 				ConvexAlignFast::FastUnrolledfwdFillMatrixLine(refSeq,qrySeq,fwdResult,readId,
+	// 					matrix->getElementUp(x - 1, y - 1)->score,
+	// 					*matrix->getElementUp(x,y - 1),
+	// 					*matrix->getElementCurr(x - 1, y),
+	// 					matrix->getElementEditCurrFast(x,y),
+	// 					*matrix->getDirectionCurrFast(x, y),
+	// 					curr_max,
+	// 					x,
+	// 					y,
+	// 					read_char_cache);
+	// 			}			
+	// 		} else {
+	// 			for(int x=xMin; x<=xFastMin; x++)
+	// 			{
+	// 				ConvexAlignFast::FastUnrolledfwdFillMatrixLine(refSeq,qrySeq,fwdResult,readId,
+	// 					matrix->getElementUp(x - 1, y - 1)->score,
+	// 					*matrix->getElementUp(x,y - 1),
+	// 					*matrix->getElementCurr(x - 1, y),
+	// 					matrix->getElementEditCurrFast(x,y),
+	// 					*matrix->getDirectionCurrFast(x, y),
+	// 					curr_max,
+	// 					x,
+	// 					y,
+	// 					read_char_cache);
+	// 			}
+
+	// 			for(int x=xFastMin+1; x<xFastMax; x++)
+	// 			{
+	// 				ConvexAlignFast::FastUnrolledfwdFillMatrixLine(refSeq,qrySeq,fwdResult,readId,
+	// 					matrix->getElementUpFast(x - 1, y - 1)->score,
+	// 					*matrix->getElementUpFast(x,y - 1),
+	// 					*matrix->getElementCurrFast(x - 1, y),
+	// 					matrix->getElementEditCurrFast(x,y),
+	// 					*matrix->getDirectionCurrFast(x, y),
+	// 					curr_max,
+	// 					x,
+	// 					y,
+	// 					read_char_cache);
+	// 			}
+
+
+	// 			for(int x=xFastMax; x<xMax; x++)
+	// 			{
+	// 				ConvexAlignFast::FastUnrolledfwdFillMatrixLine(refSeq,qrySeq,fwdResult,readId,
+	// 					matrix->getElementUp(x - 1, y - 1)->score,
+	// 					*matrix->getElementUp(x,y - 1),
+	// 					*matrix->getElementCurr(x - 1, y),
+	// 					matrix->getElementEditCurrFast(x,y),
+	// 					*matrix->getDirectionCurrFast(x, y),
+	// 					curr_max,
+	// 					x,
+	// 					y,
+	// 					read_char_cache);		
+	// 			}			
+	// 		}
+
+
+	 
+	// 		/*for (int x = std::max(0, xOffset); x < xMax; ++x) {
+	// 			if(y==0 ||
+	// 			  (x-1 <= (matrix->getCorridorOffset(y-1)) ) ||
+	// 			  (x >= (matrix->getCorridorOffset(y-1)+matrix->getCorridorLength(y-1)) )  ) {
+	// 				ConvexAlignFast::FastUnrolledfwdFillMatrixLine(refSeq,qrySeq,fwdResult,readId,
+	// 					matrix->getElementUp(x - 1, y - 1)->score,
+	// 					*matrix->getElementUp(x,y - 1),
+	// 					*matrix->getElementCurr(x - 1, y),
+	// 					matrix->getElementEditCurrFast(x,y),
+	// 					*matrix->getDirectionCurrFast(x, y),
+	// 					curr_max,
+	// 					x,
+	// 					y);
+	// 			} else {
+	// 				ConvexAlignFast::FastUnrolledfwdFillMatrixLine(refSeq,qrySeq,fwdResult,readId,
+	// 					matrix->getElementUpFast(x - 1, y - 1)->score,
+	// 					*matrix->getElementUpFast(x,y - 1),
+	// 					*matrix->getElementCurrFast(x - 1, y),
+	// 					matrix->getElementEditCurrFast(x,y),
+	// 					*matrix->getDirectionCurrFast(x, y),
+	// 					curr_max,
+	// 					x,
+	// 					y);
+	// 			}
+	// 		}*/
+
+	// 	}
+	// 	fwdResult.qend = (matrix->getHeight() - fwdResult.best_read_index) - 1;
+	// 	if (matrix->getHeight() == 0) {
+	// 		fwdResult.best_read_index = fwdResult.best_ref_index = 0;
+	// 	}
+
+	// 	return curr_max;
+	// }
+
+	// inline void FastUnrolledfwdFillMatrixLine(char const * const refSeq,
+	// 		char const * const qrySeq, FwdResults & fwdResult, int readId,
+	// 		AlignmentMatrixFast::Score diag_score,
+	// 		AlignmentMatrixFast::MatrixElement const & up,
+	// 		AlignmentMatrixFast::MatrixElement const & left,
+	// 		AlignmentMatrixFast::MatrixElement * current,
+	// 		char & currentDirection,
+	// 		AlignmentMatrixFast::Score & curr_max,
+	// 		int x,
+	// 		int y,
+	// 		char read_char_cache) {
+			
+	// 		bool const eq = read_char_cache == refSeq[x];
+	// 		AlignmentMatrixFast::Score const diag_cell = diag_score + (eq ? mat : mis);
+
+	// 		AlignmentMatrixFast::Score up_cell = 0;
+	// 		AlignmentMatrixFast::Score left_cell = 0;
+
+	// 		int ins_run = 0;
+	// 		int del_run = 0;
+
+	// 		if (up.direction == CIGAR_I) {
+	// 			ins_run = up.indelRun;
+	// 			if (up.score == 0) {
+	// 				up_cell = 0;
+	// 			} else {
+	// 				up_cell = up.score
+	// 						+ std::min(gap_ext_min,
+	// 								gap_ext + ins_run * gap_decay);
+	// 			}
+	// 		} else {
+	// 			up_cell = up.score + gap_open_read;
+	// 		}
+
+	// 		if (left.direction == CIGAR_D) {
+	// 			del_run = left.indelRun;
+	// 			if (left.score == 0) {
+	// 				left_cell = 0;
+	// 			} else {
+	// 				left_cell = left.score
+	// 						+ std::min(gap_ext_min,
+	// 								gap_ext + del_run * gap_decay);
+	// 			}
+	// 		} else {
+	// 			left_cell = left.score + gap_open_ref;
+	// 		}
+
+	// 		//find max
+	// 		AlignmentMatrixFast::Score max_cell = 0;
+	// 		max_cell = std::max(left_cell, max_cell);
+	// 		max_cell = std::max(diag_cell, max_cell);
+	// 		max_cell = std::max(up_cell, max_cell);
+
+	// 		/*int a=del_run > 0 && max_cell == left_cell;
+	// 		int b=(!a) && ins_run > 0 && max_cell == up_cell;
+	// 		int c=(!a && !b) && max_cell == diag_cell;
+	// 		int d=(!a && !b && !c) && max_cell == left_cell;
+	// 		int e=(!a && !b && !c && !d) && max_cell == up_cell;
+	// 		int f=!(a*b*c*d*e);
+
+	// 		current->score = a*max_cell;
+	// 		current->direction = a*CIGAR_D;
+	// 		currentDirection = a*CIGAR_D;
+	// 		current->indelRun = a*(del_run + 1);
+
+	// 		current->score = (!b)*current->score + b*max_cell;
+	// 		current->direction = (!b)*current->direction + b*CIGAR_I;
+	// 		currentDirection = (!b)*currentDirection + b*CIGAR_I;
+	// 		current->indelRun = (!b)*current->indelRun + b*(ins_run + 1);
+
+	// 		current->score = (!c)*current->score + c*max_cell;
+	// 		current->direction = (!c)*current->direction + (c*eq*CIGAR_EQ) + (c*(!eq)*CIGAR_X);
+	// 		currentDirection = (!c)*currentDirection + (c*eq*CIGAR_EQ) + (c*(!eq)*CIGAR_X);
+	// 		current->indelRun = (!c)*current->indelRun + c*0;
+
+	// 		current->score = (!d)*current->score + d*max_cell;
+	// 		current->direction = (!d)*current->direction + d*CIGAR_D;
+	// 		currentDirection = (!d)*currentDirection + d*CIGAR_D;
+	// 		current->indelRun = (!d)*current->indelRun + d*1;
+
+	// 		current->score = (!e)*current->score + e*max_cell;
+	// 		current->direction = (!e)*current->direction + e*CIGAR_I;
+	// 		currentDirection = (!e)*currentDirection + e*CIGAR_I;
+	// 		current->indelRun = (!e)*current->indelRun + e*1;
+
+	// 		current->score = (!f)*current->score + 0;
+	// 		current->direction = (!f)*current->direction + f*CIGAR_STOP;
+	// 		currentDirection = (!f)*currentDirection + f*CIGAR_STOP;
+	// 		current->indelRun = (!f)*current->indelRun + 0;*/
+
+	// 		if (del_run > 0 && max_cell == left_cell) {
+	// 			current->score = max_cell;
+	// 			current->direction = CIGAR_D;
+	// 			currentDirection = CIGAR_D;
+	// 			current->indelRun = del_run + 1;
+	// 		} else if (ins_run > 0 && max_cell == up_cell) {
+	// 			current->score = max_cell;
+	// 			current->direction = CIGAR_I;
+	// 			currentDirection = CIGAR_I;
+	// 			current->indelRun = ins_run + 1;
+	// 		} else if (max_cell == diag_cell) {
+	// 			current->score = max_cell;
+	// 			if (eq) {
+	// 				current->direction = CIGAR_EQ;
+	// 				currentDirection = CIGAR_EQ;
+	// 			} else {
+	// 				current->direction = CIGAR_X;
+	// 				currentDirection = CIGAR_X;
+	// 			}
+	// 			current->indelRun = 0;
+	// 		} else if (max_cell == left_cell) {
+	// 			current->score = max_cell;
+	// 			current->direction = CIGAR_D;
+	// 			currentDirection = CIGAR_D;
+	// 			current->indelRun = 1;
+	// 		} else if (max_cell == up_cell) {
+	// 			current->score = max_cell;
+	// 			current->direction = CIGAR_I;
+	// 			currentDirection = CIGAR_I;
+	// 			current->indelRun = 1;
+	// 		} else {
+	// 			current->score = 0;
+	// 			current->direction = CIGAR_STOP;
+	// 			currentDirection = CIGAR_STOP;
+	// 			current->indelRun = 0;
+	// 		}
+
+	// 		if (max_cell > curr_max) {
+	// 			curr_max = max_cell;
+	// 			fwdResult.best_ref_index = x;
+	// 			fwdResult.best_read_index = y;
+	// 			fwdResult.max_score = curr_max;
+	// 		}
+	// }
+
 	/**
 	 * Follows direction matrix (backtracking)
 	 * Creates reversed binary representation of the cigar string