Commit b638c912 authored by thomas dettbarn's avatar thomas dettbarn
Browse files

more correct calculations in right_64().

parent 81d51e9c

lib_dec/ivas_svd_dec_fx.c

+36 −13
Original line number Diff line number Diff line
@@ -1130,19 +1130,21 @@ static void biDiagonalReductionRight_64( 
    IF( LT_16( currChannel, nChannelsL ) && NE_16( currChannel, sub( nChannelsC, 1 ) ) ) /* i <=m && i !=n */

    {

        Word64 norm_64;

        Word64 abs_x;

        Word64 abs_64;

        Word16 idx;

        Word32 tmp;

        Word16 tmpe;

        Word16 iCh,jCh;

        Word32 norm_x;

        Word16 norm_x_e;

	Word32 abs_x;

	Word16 abs_x_e;

     



        idx=add(currChannel,1);



        norm_64=0;

        abs_x=0;

        abs_64=0;

        move64();

        move64();

        tmpe=add(singularVectors_e,1);
@@ -1151,12 +1153,15 @@ static void biDiagonalReductionRight_64( 
        {

            tmp=W_extract_l(W_shr(singularVectors_Left_64[currChannel][jCh],tmpe));

            norm_64=W_add(norm_64,W_mult0_32_32(tmp,tmp));

            abs_x=W_add(abs_x, W_abs(singularVectors_Left_64[currChannel][jCh]));

            abs_64=W_add(abs_64, W_abs(singularVectors_Left_64[currChannel][jCh]));

        }

        tmpe=W_norm(norm_64);

        norm_x=W_extract_h(W_shl(norm_64, tmpe ));

        norm_x_e = add(sub(shl(singularVectors_e, 1), tmpe), 3 );

	printf("%016llX --> %08X<%2X   abs_x:%016llX \x1b[0m\n",norm_64,norm_x,norm_x_e,abs_x);

        tmpe=W_abs(abs_64);

        abs_x=W_extract_h(W_shl(abs_64, tmpe ));

        abs_x_e = add(sub(shl(singularVectors_e, 1), tmpe), 3 );

	printf("%016llX --> %08X<%2X   abs_x:%016llX --> %08X<%2X\x1b[0m\n",norm_64,norm_x,norm_x_e,abs_64,abs_x,abs_x_e);

        IF (norm_x)

        {

            Word16 invVal_e;
@@ -1185,44 +1190,53 @@ static void biDiagonalReductionRight_64( 
            *g_e = L_temp_e;

            move16();

    

            tmp=W_extract_l(W_shr(singularVectors_Left_64[currChannel][idx],singularVectors_e) );

            tmp=W_extract_l(W_shr(singularVectors_Left_64[currChannel][idx],singularVectors_e+1) );

            tmpe=sub(sub(*g_e, singularVectors_e),1);

            tmpmul=W_mult0_32_32(*g,tmp);

            tmpmul=W_shl(tmpmul,tmpe);

            r_64=W_sub(tmpmul, norm_64 );

            r_64=W_sub(W_shr(tmpmul,1), W_shr(norm_64,1) );

		printf("MUL: \x1b[1;33ms:%016llX -> tmp:%08X -> g:%08X tmpmul:%016llX -> norm_64:%016llX --> R:%016llX\x1b[0m\n",singularVectors_Left_64[currChannel][idx],tmp,*g,tmpmul,norm_64,r_64);

            

            tmpe2=W_norm(r_64);

            r=W_extract_h(W_shl(r_64,tmpe2));

            r_e = sub(sub(add(add(1,singularVectors_e),*g_e),tmpe),tmpe2);



		printf("R: \x1b[1;33m%016llX --> %08X<%2X   g:%08X<%2X\x1b[0m\n",r_64,r,r_e,*g,*g_e);

    

            invVal_e = r_e;

            invVal = BASOP_Util_Divide3232_Scale_newton( MAXVAL_WORD32, r, &invVal_e);

    

            tmpe=add(32,sub(singularVectors_e,g_e));	// TODO: maybe the other way around??

            tmp64=W_shr(W_deposit32_h(g),tmpe);

//            tmpe=add(31,sub(singularVectors_e,g_e));	// TODO: maybe the other way around??

		tmpe=sub(30,singularVectors_e);

            tmp64=W_shr(W_deposit32_h(*g),tmpe);

		printf("SING: \x1b[1;33m%016llX-g:%08X (%016llX) = ",singularVectors_Left_64[currChannel][idx],*g,tmp64);

    	    singularVectors_Left_64[currChannel][idx]=W_sub(singularVectors_Left_64[currChannel][idx],tmp64);	// exponent +1

		printf("%016llX\x1b[0m\n",singularVectors_Left_64[currChannel][idx]);





            tmpe=add(singularVectors_e,2);

            tmpe=add(singularVectors_e,8);

            tmp=W_extract_l(W_shr(singularVectors_Left_64[currChannel][currChannel],tmpe));

            for ( iCh=currChannel+1; iCh<nChannelsC; iCh++)

            {

                Word32 tmp2;

                Word32 tmp3;

    

        

       		printf("norm: \x1b[1;33m"); 

                tmp2=W_extract_l(W_shr(singularVectors_Left_64[currChannel][iCh],tmpe));

                norm_64=W_mult0_32_32(tmp,tmp2);

                for (jCh=currChannel+1;jCh<nChannelsL; jCh++)

                {

                    tmp2=W_extract_l(W_shr(singularVectors_Left_64[jCh][currChannel],tmpe));

                    tmp3=W_extract_l(W_shr(singularVectors_Left_64[jCh][iCh],tmpe));

			printf("%016llx ",norm_64);

                    tmp2=W_extract_l(W_shr(singularVectors_Left_64[currChannel][jCh],tmpe));

                    tmp3=W_extract_l(W_shr(singularVectors_Left_64[iCh][jCh],tmpe));

                    norm_64=W_add(norm_64, W_mult0_32_32(tmp2,tmp3));

                }

		printf("%016llx\x1b[0m\n",norm_64);

                norm_x_e=W_norm(norm_64);

                norm_x=W_extract_h(W_shl(norm_64,norm_x_e));

                f=Mpy_32_32(norm_x,invVal);

                f_e=add(invVal_e, sub(norm_x_e, r_e ));

		printf("F: \x1b[1;33mnorm:%016llX --> %08X<%2X   %08X<%2X\x1b[0m\n",norm_64,norm_x,norm_x_e,f,f_e);

    

                for (jCh=idx;jCh<nChannelsL; jCh++)

                {
@@ -1593,7 +1607,7 @@ static void biDiagonalReductionRight_fx( 
            norm_x = BASOP_Util_Add_Mant32Exp( norm_x, norm_x_e, Mpy_32_32( singularVectors[currChannel][jCh], singularVectors[currChannel][jCh] ), shl( singularVectors2_e[currChannel][jCh], 1 ), &norm_x_e ); /* exp(norm_x_e) */

            abs_x = BASOP_Util_Add_Mant32Exp( abs_x, abs_x_e, L_abs( singularVectors[currChannel][jCh] ), singularVectors2_e[currChannel][jCh], &abs_x_e );

        }

	printf("%08X<%2X   %08X<%2X \x1b[0m\n",norm_x,norm_x_e,abs_x,abs_x_e);

	printf("%08X<%2X   abs_x:%08X<%2X \x1b[0m\n",norm_x,norm_x_e,abs_x,abs_x_e);



        IF( norm_x ) /*(fabsf(*sig_x) > EPSILON * fabsf(*sig_x)) { */

        {
@@ -1618,9 +1632,13 @@ static void biDiagonalReductionRight_fx( 




            r = BASOP_Util_Add_Mant32Exp( Mpy_32_32( ( *g ), singularVectors[currChannel][idx] ), add( *g_e, singularVectors2_e[currChannel][idx] ), -norm_x, norm_x_e, &r_e );                                      /* exp(r_e) */

		printf("MUL: s:%08X -> g:%08X ->tmpmul:%08X norm:%08X --> r:%08X\n",singularVectors[currChannel][idx],*g,Mpy_32_32( ( *g ), singularVectors[currChannel][idx]),norm_x,r);

		printf("SING: %08X-%08X = ",singularVectors[currChannel][idx],*g);

            singularVectors[currChannel][idx] = BASOP_Util_Add_Mant32Exp( singularVectors[currChannel][idx], singularVectors2_e[currChannel][idx], -( *g ), ( *g_e ), &singularVectors2_e[currChannel][idx] ); /* exp(sing_exp) */

		printf("%08X\n",singularVectors[currChannel][idx]);

            move32();



		printf("R: --> %08X<%2X g:%08X<%2x\n",r,r_e,*g,*g_e);



            FOR( iCh = currChannel + 1; iCh < nChannelsL; iCh++ ) /*  nChannelsL */

            {
@@ -1628,12 +1646,17 @@ static void biDiagonalReductionRight_fx( 
                move32();

                norm_x_e = 0;

                move16();

		printf("norm: ");

                FOR( jCh = idx; jCh < nChannelsC; jCh++ ) /* nChannelsC */

                {

                    norm_x = BASOP_Util_Add_Mant32Exp( norm_x, norm_x_e, Mpy_32_32( singularVectors[iCh][jCh], singularVectors[currChannel][jCh] ), add( singularVectors2_e[iCh][jCh], singularVectors2_e[currChannel][jCh] ), &norm_x_e ); /* exp(norm_x_e) */

			printf("%08X<%02X  ",norm_x,norm_x_e);

                }

		printf("\n");

		printf("F: norm %08X<%2X  ",norm_x,norm_x_e);

                norm_x = BASOP_Util_Divide3232_Scale_newton( norm_x, r, &invVal_e );

                norm_x_e = add( invVal_e, sub(norm_x_e, r_e ) );

		printf("-> %08X<%2X\n",norm_x,norm_x_e);

                FOR( jCh = idx; jCh < nChannelsC; jCh++ ) /*  nChannelsC */

                {

                    singularVectors[iCh][jCh] = BASOP_Util_Add_Mant32Exp( singularVectors[iCh][jCh], singularVectors2_e[iCh][jCh], Mpy_32_32( norm_x, singularVectors[currChannel][jCh] ), add( norm_x_e, singularVectors2_e[currChannel][jCh] ), &singularVectors2_e[iCh][jCh] ); /* exp(sing_exp2) */