format pe gui 4.0
include 'win32ax.inc'

section '' code executable import readable writeable
  library kernel32,'kernel32.dll',\
          user32,'user32.dll'
  include 'api\kernel32.inc'
  include 'api\user32.inc'






;tf f1,$3fff-0000,$8000000000000000 ;2^0*1
;tf f2,$bfff-0064,$8000000000000000 ;2^-64*1
macro tf name*, exponent*, significand*{
  label name tbyte at $
  dq significand
  dw exponent
}

;include '../macro.inc'
;include '../proc.inc'



;buf80bytes db 80 dup 0
f dt -65535.0
;tf f,$403e,$ffffffffffffffff


entry $
        mov     eax,f
        call    f2xey
        ret



/* IT MUST BE TESTED
  Prepares floating value for text formating.  This is an intermediate procedure.
  eax         -> tbyte float location
  ecx:ebx:eax <- decimal mantissa given as a float, absolute value
  edx         <- bit 31..31: sign of decimal mantissa, 1='-', 0='+'
                 bit 30..30: sign of decimal exponent, 1='-', 0='+'
                 bit 00..15: decimal exponent given as an integer, absolute value
  flags.cf    <- 0 if not zero
                 1 if zero, ecx.31 equals zero's sign, ebx:eax=0, edx=?
  flags       <- ?
*/
proc f2xey
        ;
        ;decom`pose number
        ;
        stc                          ;ieee stupido mode ON
        call    fxtract_
        jnc     .not_zero
        ret     0

.not_zero:
        ;
        ;10^x = 2^unbisased * normalized      =
        ;     = 2^unb * 2^log(2;norm)         =
        ;     = 2^(unb + log(2;norm))        ->
        ;   x = log(10;2^(unb + log(2;norm))) =
        ;     = log(10;2)*(unb + log(2;norm))
        ;
        ;log(10;2)       = 3FFD 9A209A84 FBCFF799
        ;                =     .4D104D42'7DE7FBCC'8xxxxxxx
        ;unb             = -16'445..+16'384
        ;log(2;norm)     = [0..1)
        ;min log(2;norm) = log(2;$3fff'80000000'00000001)
        ;                = 3FC0 B8AA3B29 5C17F0BB ->
        ;               -> ie $3fff-$3fc0=$3f=63
        ;Performing "unb + log(2;norm)", fpu loses a lot of fractional bits.
        ;That's why we are to use 128 bits for fraction.
        ;Mention: for .fraction only, hence 3FFD 9A209A84 FBCFF799 =
        ;=9A209A84 FBCFF799 shr 1 but not 2
        ;
        push    esi edi ebp

        ;
        ;get log(2;normalized) into ST0
        ;
        fld1
        push    word $3fff ebx eax
        fld     tbyte[esp]
        fyl2x                        ;got = [0..1)

        ;
        ;unnormalize ST0 into BP:.EBX:EAX:ESI:EDI as common whole.fraction
        ;
        fstp    tbyte[esp]
        mov     ebp,ecx              ;keep value sign(ebp.31) and unb(ebp.00..15)
        pop     eax ebx cx           ;operate with registers
        neg     cx                   ;getting shift count to get actual value
        jnc     .idle                ;log(2;norm)=0
        add     cx,$3fff-1           ;-1..-62
.idle:  xor     esi,esi              ;complete fraction
        xor     edi,edi              ;complete fraction
        cmp     cl,32
        jb      .b32                 ;small shift
        xchg    eax,esi
        xchg    ebx,eax
        sub     cl,32                ;denormalized partially
.b32:   shrd    edi,esi,cl
        shrd    esi,eax,cl
        shrd    eax,ebx,cl
        shr     ebx,cl               ;denormalized fully - actual log(2; normalized) got

        ;
        ;add unbiased in BP to log(2;normalized) in .EBX:EAX:ESI:EDI
        ;
        test    bp,bp
        jns     .unbiased_added      ;whole.fraction :)^
        neg     bp                   ;how would you do with -1+0.3? -> -(1-0.3)=-0.7
        bts     ebp,30               ;as floats stored positive, they need a bit saying "-" or "+"
        neg     ebx                  ;subtract
        sbb     ebp,0                ;and chain; we could also use "sbb bp,0", nah
        neg     eax
        sbb     ebx,0
        sbb     ebp,0
        neg     esi
        sbb     eax,0
        sbb     ebx,0
        sbb     ebp,0
        neg     edi
        sbb     esi,0
        sbb     eax,0
        sbb     ebx,0
        sbb     ebp,0                ;negative whole.fraction got
.unbiased_added:


        ;
        ;multiply BP:.EBX:EAX:ESI:EDI by log(10;2)
        ;
        push    ebx eax esi edi
        mov     esi,esp              ;<_128mul128>.m2
        mov     ecx,$4D10'4D42
        mov     ebx,$7DE7'FBCC
        mov     eax,$8000'0000
        push    ecx ebx eax 0
        mov     edi,esp              ;<_128mul128>.m1
        stdcall _128mul128,edi,esi
        movzx   esi,bp               ;we've multiplied .fractions, now whole.*log(10;2)
        xor     bp,bp                ;clean whole. part
        mul     esi
        add     [edi+5*4],eax
        adc     [edi+6*4],edx
        adc     dword[edi+7*4],0
        adc     ebp,0
        mov     eax,ebx
        mul     esi
        add     [edi+6*4],eax
        adc     [edi+7*4],edx
        adc     ebp,0
        mov     eax,ecx
        mul     esi
        add     [edi+7*4],eax
        adc     ebp,edx              ;got whole.fraction

        ;
        ;getting 10^.fraction = (10^-1..10^+1) = (0.1..10)
        ;
        std                          ;searching for the highest non-zero dword...
        add     edi,7*4              ;assume we've found it
        xor     eax,eax
        mov     ecx,8                ;256bits=8dwords
        repe    scasd
        cld
        jne     .fraction_nz         ;found, else .fraction zero
                                     ;10^0=1...
        mov     ecx,$3fff            ;exponent=0 ie 2^0
        mov     ebx,$8000'0000       ;normalized 1.0
        xor     eax,eax              ;too
        mov     edx,ebp              ;edx.31=input sign, edx.30=10's integral power sign, edx.00.15=10's integral power
        add     esp,8*4              ;get rid of 256bit .fraction
        jmp     .tired_edemko        ;bye









































.fraction_nz:
hlt
        mov     ch,cl
        xor     esi,esi
        mov     ebx,[edi+4]
        test    cl,cl
        jz      .normalize_fraction
        mov     eax,[edi]
        dec     cl
        jz      .normalize_fraction
        mov     esi,[edi-4]
.normalize_fraction:
        neg     ch
        add     ch,7
        shl     ch,5
        bsr     edx,ebx
        neg     edx
        add     edx,31
        add     ch,dl
        mov     cl,dl
        shld    ebx,eax,cl
        shld    eax,esi,cl
        movzx   ecx,ch
        neg     ecx
        add     ecx,$3fff-1
        shl     ecx,1
        bt      ebp,30
        rcr     cx,1
        push    ecx ebx eax
        fld     tbyte[esp]
        call    power10
hlt



























.tired_edemko:
        pop     ebp edi esi
        clc
        ret     0
endp

























































;/*
;  address of tbyte float -> eax
;  flags.cf               -> ignored if qword[eax]=0
;                            else it defines how to process a $0000'xxxxxxxx'xxxxxxxx case
;                            if cf=0: no stupido IEEE processing will be done
;                                     this is not recommended if no exotic wanted
;                                     $0000'8xxxxxxx'xxxxxxxx ->fxtract-> cx = -16'383
;                                     $0000'4xxxxxxx'xxxxxxxx ->fxtract-> cx = -16'384
;                            if cf=1: increase exponent by 1 like fpu's FXTRACT does
;                                     $0000'8xxxxxxx'xxxxxxxx ->fxtract-> cx = -16'382
;                                     $0000'4xxxxxxx'xxxxxxxx ->fxtract-> cx = -16'383
;  ebx:eax                <- normalized mantissa
;  ecx                    <- bit 31   : sign of mantissa
;                            bit 0..15: unbiased integer exponent within -16'446(MINUS stupido.ieee.flags.cf)..+16'384
;                            bit rest : 0
;  flags.cf               <- 1 if mantissa was zero
;  ------------------------------------------------------
;  possible usage:
;    2^unbiased * normalized = 10^x
;    x = log(10;2^unbiased * normalized)
;      = log(10;2^unbiased) + log(10;normalized)
;      = log(10;2)*unbiased + log(10;2)*log(2;normalized)
;      = log(10;2)*(unbiased + log(2;normalized))
;      = Frank Airline in my dreams, thank you
;  possible usage:
;    2^unbiased = 10^x
;    x=log(10;2^unbiased)=log(10;2)*unbiased
;  possible usage:
;    10^x = 2^(unbiased + log(2;normalized))
;    x    = log(10;2) * (unbiased + log(2;normalized))
;  possible usage:
;    boundless division and multiplication
;*/
proc fxtract_
        push    edx
        pushf                    ;save ieee.stupido
        movzx   ecx,word[eax+8]
        shl     cx,1             ;drop sign bit
        rcr     ecx,1            ;into the top
        sub     cx,$3fff         ;to easily get unbiased exponent
        mov     ebx,[eax+4]      ;mantissa
        mov     eax,[eax]        ;mantissa
        mov     edx,eax
        or      edx,ebx
        jnz     .nz              ;else zero found
        or      byte[esp],1      ;flags.cf=1
.done:  popf
        pop     edx
        ret     0                ;fasmme
.nz:    btr     dword[esp],0     ;flags.cf=0, get ieee.stupido
        jnc     .no_ieee         ;no ieee.stupido
        cmp     cx,$c001         ;have we met a stupido condition(-16'383)?
        jne     .no_ieee         ;no
        inc     ecx              ;do stupido
        test    ebx,ebx
        js      .done            ;there is nothing to normalize in mantissa, return -16'383
.no_ieee:
        test    ebx,ebx
        jnz     .ebx_nz          ;end normalization
        xchg    ebx,eax
        sub     cx,32            ;partial normalization done
.ebx_nz:bsr     edx,ebx          ;get times to normalize
        neg     edx
        add     edx,31
        sub     cx,dx            ;_balance_
        xchg    cl,dl
        shld    ebx,eax,cl
        shl     eax,cl
        mov     cl,dl
        jmp     .done            ;2010_08_21
endp




;/*
;  [m2:m1] = [m2]*[m1].
;  flags   <-?
;  low 128 bits stored at <m1> then high 128 bits stored at <m2>
;  there will be fun when buffers overlap :)
;*/
proc _128mul128; m1,m2
        push    ebp
        mov     ebp,esp
        push    edi
        mov     edi,[ebp+2*4]           ;m1
        push    esi
        mov     esi,[ebp+3*4]           ;m2
        push    edx ecx ebx eax         ;we'll use stepping _64mul64
                                        ;[m1]*[m2] =  [m1.1]*[m2.1]          + ;00I
                                        ;          + ([m1.1]*[m2.2]) shl 064 + ;0II
                                        ;          + ([m1.2]*[m2.1]) shl 064 + ;III
                                        ;          + ([m1.2]*[m2.2]) shl 128   ;0IV
        mov     edx,[esi+4]
        mov     ecx,[esi]
        mov     ebx,[edi+4]
        mov     eax,[edi]
        call    _64mul64
        push    0 0 0 0 edx ecx ebx eax ;00I
        mov     edx,[esi+12]
        mov     ecx,[esi+8]
        mov     ebx,[edi+4]
        mov     eax,[edi]
        call    _64mul64
        add     [ebp-12*4],eax
        adc     [ebp-11*4],ebx
        adc     [ebp-10*4],ecx
        adc     [ebp-9*4],edx
        adc     dword[ebp-8*4],0        ;0II
        mov     edx,[esi+4]
        mov     ecx,[esi]
        mov     ebx,[edi+12]
        mov     eax,[edi+8]
        call    _64mul64
        add     [ebp-12*4],eax
        adc     [ebp-11*4],ebx
        adc     [ebp-10*4],ecx
        adc     [ebp-9*4],edx
        adc     dword[ebp-8*4],0        ;III
        mov     edx,[esi+12]
        mov     ecx,[esi+8]
        mov     ebx,[edi+12]
        mov     eax,[edi+8]
        call    _64mul64
        add     [ebp-10*4],eax
        adc     [ebp-9*4],ebx
        adc     [ebp-8*4],ecx
        adc     [ebp-7*4],edx           ;0IV
        pop     dword[edi] dword[edi+4] dword[edi+8] dword[edi+12] dword[esi] dword[esi+4] dword[esi+8] dword[esi+12]

        pop     eax ebx ecx edx
        pop     esi edi
        leave
        ret     8                       ;2010_08_20

endp




;you may skip it though it's interesting
; $ffffffff(edx)ffffffff(ecx) * $ffffffff(ebx)ffffffff(eax)      =
;=$ffffffff'ffffffff * $1'00000000'00000000 - $ffffffff'ffffffff =
;=$ffffffff'ffffffff'00000000'00000000      - $ffffffff'ffffffff =
;=$ffffffff'ffffffff'00000000'00000000  - 1 - $ffffffff'fffffffe =
;=$ffffffff'fffffffe'ffffffff'ffffffff      - $ffffffff'fffffffe =
;=$ffffffff'fffffffe'00000000'00000001, 128 bits :b

;edx:ecx:ebx:eax = edx:ecx*ebx:eax = ecx*ebx:eax+(edx*ebx:eax)shl 32
;flags           = ?
proc _64mul64
        push    eax ebx edx
        call    _64mul32
        xchg    eax,[esp+8]
        xchg    ebx,[esp+4]
        xchg    ecx,[esp]
        call    _64mul32
        mov     edx,ecx
        mov     ecx,ebx
        mov     ebx,eax
        pop     eax
        add     ecx,eax
        adc     edx,0
        pop     eax
        add     ebx,eax
        adc     ecx,0
        adc     edx,0
        pop     eax
        ret     0
endp




;ecx:ebx:eax = ecx*ebx:eax = ecx*ebx:0+ecx*0:eax
;edx         = ?, ecx*0:eax shr 32 actually
;flags       = ?
proc _64mul32
        xchg    eax,ebx
        mul     ecx
        xchg    ebx,eax
        xchg    ecx,edx
        mul     edx
        add     ebx,edx
        adc     ecx,0
        ret     0
endp




;power   -> st0
;st0     <- 10^power
;fpu.tos <- returned
;ensure st5..st7 are free before calling this proc
proc power10
                       ;10^power=2^x
                       ;x=log(2;10^power)=log(2;10)*power
        fldl2t         ;log(2;10)
        fmulp          ;x

                       ;2^x=2^x.fraction*2^x.whole=2^x.fraction shl x.whole
        fld     st0    ;x
        frndint        ;x.whole
        fsub    st1,st0;x.fraction
        fxch
        f2xm1          ;2^x.fraction -1
        fld1
        faddp          ;2^x.fraction
        fscale         ;shl x.whole
        fstp    st1    ;10^x
        ret     0      ;-)
endp
    

fld [t0_256]
fldl2t
fmulp
f2xm1 ; this is okay because -1 < 0.256 * log2(10) < 1    

;/*
;  Recognize hex-string into a tbyte float.
;  src         = string
;  options     = bits 00..007: dummy char not causing errors
;  flags.cf    <-error
;  ecx:ebx:eax <- result
;*/
proc bsa_h2f; src, options
        stc
        pushf                    ;error=1
        push    esi
        mov     esi,[esp+3*4]    ;src
        test    esi,esi
        jz      .done            ;src=0
        movzx   ebx,byte[esi-4]
        lea     esi,[esi+ebx-1]  ;switch LeastSignificantByte
        shl     ebx,8
        jz      .done            ;src.body_len=0 | src.body_len>255
        mov     bl,[esp+4*4]     ;options(dummy)
        xor     ecx,ecx
        push    ecx ecx ecx      ;result
        std
.load:  lodsb
        cmp     al,bl
        je      .continue        ;dummy
        call    ASCII
        sub     al,'0'
        jz      .register        ;0*smth=0
        cmp     ecx,80-4
        ja      .pop             ;over(80bits)flow
        cmp     al,10
        jb      .input           ;00..10
        cmp     al,'F'-'0'
        ja      .pop             ;not hex
        cmp     al,'A'-'0'
        jb      .pop             ;not hex
        sub     al,7             ;10..15
.input: movzx   eax,al
        cmp     cl,64
        jb      .32?             ;we are in 64..79
        sub     cl,64
        shl     eax,cl
        add     cl,64
        or      [esp+2*4],eax    ;update result.64..79
        jmp     .register
.32?:   cmp     cl,32
        jb      .32              ;we are in 32..63
        sub     cl,32
        shl     eax,cl
        add     cl,32
        or      [esp+1*4],eax    ;update result.32..63, .., i've got used counting from 1
        jmp     .register
.32:    shl     eax,cl           ;we are in 0..31
        or      [esp],eax
.register:
        add     ecx,4
.continue:
        dec     bh
        jnz     .load            ;till the end of the sun
        test    ecx,ecx
        jz      .pop             ;dummy string
        dec     byte[esp+(3+1)*4];drop flags.cf
.pop:   pop     eax ebx ecx
.done:  pop     esi
        popf
        ret     2*4              ;2010_08_09
endp
    

   10^.fraction_256                                =
 = wait, let us limit it to .fraction128           =
 = 10^.fraction[127..064] * 10^.fraction[063..000] =
 = happy                                           =
 = multiply by 10^19 or 10^20 to get 19 digits of fraction

 10^19 and 10^20, suppose those are 10^1 and 10^2 ->
 1.0 * 10^1 = 10
 0.1 * 10^2 = 10
 =thank you
    

;
;getting power of ten we use a formula: 10^x=2^y, y=log(2;10)*x
;finally when fpu executes F2XM1 we must add 1 where rounding occurs
;hence a new formula, assuming that .fraction takes 128 bits:
;      1_F2XM1*2^1_FSCALE * 2_F2XM1*2^2_FSCALE         =
;->    removing M1 factor for very member              ->
;->    1_F2XM1*2^(1_FSCALE+1) * 2_F2XM1*2^(2_FSCALE+1) =
;=     1_F2XM1*2_F2XM1*2^(1_FSCALE+2_FSCALE+2)         eg
;eg    (2-1)*2^0 * (2-1)*2^0 =
;=     1*2^0     * 1*2^0     ->
;->    1*1*2^(0+0+2)         =
;=     4                     check
;check 2*2^0 * 2*2^0         =
;=     4                     =
;=     :)
;