Colorize module

wilbert · Post by **wilbert** » Thu Feb 04, 2016 8:00 am

A SSE2 optimized module to colorize an image, inspired by the Colorize an image thread from BasicallyPure.

The easiest way to use is Colorize::ColorizeImage(0, $0080c0, Colorize::#Multiply)
The effect of the Colorize::#Overlay mode is similar to what is called Duotone.

ColorizeImageSMH allows you to set custom colors for shadow, midtone and highlight.
It's slower but gives you a lot more creative options. Two examples
Colorize::ColorizeImageSMH(0, $800000, $907040, $a0e080)
Colorize::ColorizeImageSMH(0, $800000, $ff80ff, $ffffff)

It is strongly recommended to use 32 bit images since they are processed much faster as 24 bit images.

Code: Select all

; Colorize module by Wilbert

; Latest update Feb 6, 2016

; SSE2 required

; Coefficients used for luma conversion: 0.2126, 0.7152, and 0.0722 (BT.709)

DeclareModule Colorize
  
  Enumeration
    #Multiply
    #Overlay
    #Screen
  EndEnumeration
  
  Declare ColorizeImage(Image, Color.l = $ffffff, BlendMode = #Multiply)
  Declare ColorizeImageSMH(Image, Shadow = $000000, Midtone = $808080, Highlight = $ffffff)
  
  Declare BMMultiply  (*PixArrayIn, *PixArrayOut, PixCount.l , Color.l = $ffffffff, SwapRB = #False)
  Declare BMOverlay   (*PixArrayIn, *PixArrayOut, PixCount.l , Color.l = $ff808080, SwapRB = #False)
  Declare BMScreen    (*PixArrayIn, *PixArrayOut, PixCount.l , Color.l = $ff000000, SwapRB = #False)
  
  Declare SMH(*PixArrayIn, *PixArrayOut, PixCount, Shadow = $000000, Midtone = $808080, Highlight = $ffffff, SwapRB = #False)
  
EndDeclareModule

Module Colorize
  
  DisableDebugger
  EnableExplicit
  EnableASM
  
  CompilerIf #PB_Compiler_Processor = #PB_Processor_x86
    Macro rax : eax : EndMacro
    Macro rdx : edx : EndMacro
  CompilerEndIf 
  
  Procedure ColorizeImage(Image, Color.l = $ffffff, BlendMode = #Multiply)
    Protected.i x, y, w, h, *db, pf, srb
    If IsImage(Image)
      Color | $ff000000
      w = ImageWidth(Image)
      h = ImageHeight(Image)
      StartDrawing(ImageOutput(Image))
      pf = DrawingBufferPixelFormat()
      
      If pf & (#PB_PixelFormat_32Bits_RGB | #PB_PixelFormat_32Bits_BGR)
        
        ; process pixel buffer
        If pf & #PB_PixelFormat_32Bits_BGR
          srb = #True
        EndIf
        *db = DrawingBuffer()
        Select BlendMode
          Case #Overlay : BMOverlay   (*db, *db, h * DrawingBufferPitch() >> 2, Color, srb)
          Case #Screen  : BMScreen    (*db, *db, h * DrawingBufferPitch() >> 2, Color, srb)
          Default       : BMMultiply  (*db, *db, h * DrawingBufferPitch() >> 2, Color, srb)           
        EndSelect        
        
      Else
        
        DrawingMode(#PB_2DDrawing_AllChannels)
        Protected Dim PixelLine.l(w - 1) 
        While y < h
          ; read a line of pixels
          x = 0
          While x < w
            PixelLine(x) = Point(x, y)
            x + 1
          Wend
          ; process a line of pixels
          Select BlendMode
            Case #Overlay : BMOverlay   (@PixelLine(), @PixelLine(), w, Color)
            Case #Screen  : BMScreen    (@PixelLine(), @PixelLine(), w, Color)
            Default       : BMMultiply  (@PixelLine(), @PixelLine(), w, Color)           
          EndSelect
          ; write a line of pixels
          x = 0
          While x < w
            Plot(x, y, PixelLine(x))
            x + 1
          Wend
          y + 1
        Wend
        
      EndIf
      StopDrawing()
    EndIf  
  EndProcedure
  
  Procedure ColorizeImageSMH(Image, Shadow = $000000, Midtone = $808080, Highlight = $ffffff)
    Protected.i x, y, w, h, *db, pf, srb
    If IsImage(Image)
      w = ImageWidth(Image)
      h = ImageHeight(Image)
      StartDrawing(ImageOutput(Image))
      pf = DrawingBufferPixelFormat()
      
      If pf & (#PB_PixelFormat_32Bits_RGB | #PB_PixelFormat_32Bits_BGR)
        
        ; process pixel buffer
        If pf & #PB_PixelFormat_32Bits_BGR
          srb = #True
        EndIf
        *db = DrawingBuffer()
        SMH(*db, *db, h * DrawingBufferPitch() >> 2, Shadow, Midtone, Highlight, srb)
        
      Else
        
        DrawingMode(#PB_2DDrawing_AllChannels)
        Protected Dim PixelLine.l(w - 1) 
        While y < h
          ; read a line of pixels
          x = 0
          While x < w
            PixelLine(x) = Point(x, y)
            x + 1
          Wend
          ; process a line of pixels
          SMH(@PixelLine(), @PixelLine(), w, Shadow, Midtone, Highlight)
          ; write a line of pixels
          x = 0
          While x < w
            Plot(x, y, PixelLine(x))
            x + 1
          Wend
          y + 1
        Wend
        
      EndIf
      StopDrawing()
    EndIf  
  EndProcedure
  
  Procedure BMMultiply(*PixArrayIn, *PixArrayOut, PixCount.l , Color.l = $ffffffff, SwapRB = #False)
    ; load pixel array pointers and count
    mov rax, *PixArrayIn
    mov rdx, *PixArrayOut
    mov ecx, PixCount
    ; load luma constant and color
    movd xmm2, Color  
    !movq xmm1, [colorize.l_mluma_rgba]
    !punpcklbw xmm2, xmm2
    ; swap red and blue channels if requested
    bt SwapRB, 0
    !jnc colorize.l_mloop
    !pshuflw xmm1, xmm1, 11000110b
    !pshuflw xmm2, xmm2, 11000110b
    ; main loop
    !colorize.l_mloop:
    ; load pixel color
    movd xmm0, [rax]
    !punpcklbw xmm0, xmm0
    ; convert pixel color to luma
    !pmulhuw xmm0, xmm1
    !pshuflw xmm3, xmm0, 11010010b
    !paddw xmm0, xmm3
    !pshuflw xmm3, xmm3, 11010010b
    !paddw xmm0, xmm3
    ; multiply with color
    !pmulhuw xmm0, xmm2
    ; convert back to 32 bit color value
    !psrlw xmm0, 8
    !packuswb xmm0, xmm0
    movd [rdx], xmm0
    add rax, 4
    add rdx, 4
    !sub ecx, 1
    !jnz colorize.l_mloop
    ProcedureReturn
    !colorize.l_mluma_rgba: dq 0x5555127cb717366d ; aaaa-bbbb-gggg-rrrr
  EndProcedure
  
  Procedure BMOverlay(*PixArrayIn, *PixArrayOut, PixCount.l , Color.l = $ff808080, SwapRB = #False)
    ; load pixel array pointers and count
    mov rax, *PixArrayIn
    mov rdx, *PixArrayOut
    mov ecx, PixCount  
    ; load luma constant and color
    movd xmm2, Color
    !movq xmm1, [colorize.l_oluma_rgba]
    !punpcklbw xmm2, xmm2
    ; swap red and blue channels if requested
    bt SwapRB, 0
    !jnc colorize.l_oloop
    !pshuflw xmm1, xmm1, 11000110b
    !pshuflw xmm2, xmm2, 11000110b
    ; main loop
    !colorize.l_oloop:
    ; load pixel color
    movd xmm0, [rax]
    !punpcklbw xmm0, xmm0
    ; convert pixel color to luma
    !pmulhuw xmm0, xmm1
    !pshuflw xmm3, xmm0, 11010010b
    !paddw xmm0, xmm3
    !pshuflw xmm3, xmm3, 11010010b
    !paddw xmm0, xmm3
    ; create inversion mask
    !pxor xmm5, xmm5
    !pcmpgtw xmm5, xmm0
    ; invert color and luma if needed
    !movdqa xmm3, xmm2
    !pxor xmm0, xmm5
    !pxor xmm3, xmm5
    ; luma * 2
    !paddw xmm0, xmm0
    ; multiply
    !pmulhuw xmm0, xmm3
    ; invert back result if needed
    !pxor xmm0, xmm5
    ; convert back to 32 bit color value
    !psrlw xmm0, 8
    !packuswb xmm0, xmm0
    movd [rdx], xmm0
    add rax, 4
    add rdx, 4
    !sub ecx, 1
    !jnz colorize.l_oloop
    ProcedureReturn
    !colorize.l_oluma_rgba: dq 0x2aaa127cb717366d ; aaaa-bbbb-gggg-rrrr
  EndProcedure
  
  Procedure BMScreen(*PixArrayIn, *PixArrayOut, PixCount.l , Color.l = $ff000000, SwapRB = #False)
    ; load pixel array pointers and count
    mov rax, *PixArrayIn
    mov rdx, *PixArrayOut
    mov ecx, PixCount  
    ; load luma constant and color
    movd xmm2, Color
    !movq xmm1, [colorize.l_sluma_rgba]
    !punpcklbw xmm2, xmm2
    ; create inversion mask
    !pcmpeqw xmm5, xmm5
    !psrlq xmm5, 16
    ; invert color
    !pxor xmm2, xmm5
    ; swap red and blue channels if requested
    bt SwapRB, 0
    !jnc colorize.l_sloop    
    !pshuflw xmm1, xmm1, 11000110b
    !pshuflw xmm2, xmm2, 11000110b
    ; main loop
    !colorize.l_sloop:
    ; load pixel color
    movd xmm0, [rax]
    !punpcklbw xmm0, xmm0
    ; convert pixel color to luma
    !pmulhuw xmm0, xmm1
    !pshuflw xmm3, xmm0, 11010010b
    !paddw xmm0, xmm3
    !pshuflw xmm3, xmm3, 11010010b
    !paddw xmm0, xmm3
    ; invert luma
    !pxor xmm0, xmm5
    ; multiply with inverted color
    !pmulhuw xmm0, xmm2
    ; invert back to 32 bit color value
    !pxor xmm0, xmm5
    !psrlw xmm0, 8
    !packuswb xmm0, xmm0
    movd [rdx], xmm0
    add rax, 4
    add rdx, 4
    !sub ecx, 1
    !jnz colorize.l_sloop
    ProcedureReturn
    !colorize.l_sluma_rgba: dq 0x5555127cb717366d ; aaaa-bbbb-gggg-rrrr
  EndProcedure
  
  Procedure SMH(*PixArrayIn, *PixArrayOut, PixCount, Shadow = $000000, Midtone = $808080, Highlight = $ffffff, SwapRB = #False)
    CompilerIf #PB_Compiler_OS = #PB_OS_Windows And #PB_Compiler_Processor = #PB_Processor_x64
      !movups [rsp - 16], xmm6  ; backup xmm6 (required for Win64)
    CompilerEndIf    
    ; load pixel array pointers and count
    mov rax, *PixArrayIn
    mov rdx, *PixArrayOut
    !mov ecx, [p.v_PixCount]
    ; load luma constant and color
    !movq xmm1, [colorize.l_smhluma_rgba]
    !movd xmm2, [p.v_Shadow]
    !movd xmm3, [p.v_Midtone]
    !movd xmm4, [p.v_Highlight]
    !pcmpeqw xmm5, xmm5
    !psllq xmm5, 48
    !punpcklbw xmm2, xmm2
    !punpcklbw xmm3, xmm3
    !punpcklbw xmm4, xmm4
    !por xmm2, xmm5
    !por xmm3, xmm5
    !por xmm4, xmm5
    !pxor xmm4, xmm2
    !pshuflw xmm6, xmm5, 00111111b
    ; swap red and blue channels if requested
    bt SwapRB, 0
    !jnc colorize.l_smhloop
    !pshuflw xmm1, xmm1, 11000110b
    !pshuflw xmm2, xmm2, 11000110b
    !pshuflw xmm3, xmm3, 11000110b
    !pshuflw xmm4, xmm4, 11000110b
    ; main loop
    !colorize.l_smhloop:
    ; load pixel color
    movd xmm0, [rax]
    !punpcklbw xmm0, xmm0
    ; convert pixel color to luma
    !pmulhuw xmm0, xmm1
    !pshuflw xmm5, xmm0, 11010010b
    !paddw xmm0, xmm5
    !pshuflw xmm5, xmm5, 11010010b
    !paddw xmm0, xmm5
    ; create mask
    !pxor xmm5, xmm5
    !pcmpgtw xmm5, xmm0
    ; change luma range
    !paddw xmm0, xmm0
    !pxor xmm0, xmm5
    !pxor xmm0, xmm6
    ; select shadow or highlight color
    !pand xmm5, xmm4
    !pxor xmm5, xmm2
    ; blend with mid color
    !pmulhuw xmm5, xmm0
    !pxor xmm0, xmm6
    !pmulhuw xmm0, xmm3
    !paddw xmm0, xmm5
    ; convert back to 32 bit color value
    !psrlw xmm0, 8
    !packuswb xmm0, xmm0
    movd [rdx], xmm0
    add rax, 4
    add rdx, 4
    !sub ecx, 1
    !jnz colorize.l_smhloop
    CompilerIf #PB_Compiler_OS = #PB_OS_Windows And #PB_Compiler_Processor = #PB_Processor_x64
      !movups xmm6, [rsp - 16]
    CompilerEndIf    
    ProcedureReturn
    !colorize.l_smhluma_rgba: dq 0x1555127cb717366d ; aaaa-bbbb-gggg-rrrr
  EndProcedure  
  
EndModule

walbus · Post by **walbus** » Thu Feb 04, 2016 8:51 am

Wow - many thanks Wildbret !!!

wilbert · Post by **wilbert** » Fri Feb 05, 2016 4:02 pm

I added a ColorizeImageSMH procedure.
It allows you to set separate colors for shadow, midtone and highlight.

For those who don't want to use asm or just want to understand the ColorizeSMH procedure, here it is in PureBasic code only

Code: Select all

Procedure ColorizeImageSMH(Image, Shadow = $000000, Midtone = $808080, Highlight = $ffffff)
  
  Protected.i c, x, y, w, h
  Protected.u l0,l1, r0,g0,b0, r1,g1,b1, r2,g2,b2
  
  ; expand color components to 16 bit unsigned
  r0 = Red  (Shadow) * 257
  g0 = Green(Shadow) * 257
  b0 = Blue (Shadow) * 257
  r1 = Red  (Midtone) * 257
  g1 = Green(Midtone) * 257
  b1 = Blue (Midtone) * 257
  r2 = Red  (Highlight) * 257
  g2 = Green(Highlight) * 257
  b2 = Blue (Highlight) * 257
  
  ; process image
  If IsImage(Image) And StartDrawing(ImageOutput(Image))
    DrawingMode(#PB_2DDrawing_AllChannels)
    w = OutputWidth()
    h = OutputHeight()
    While y < h
      x = 0
      While x < w
        c = Point(x, y)
        ; calculate luma
        l0 = ((c & $ff)*$36a361 + (c >> 8 & $ff)*$b7ce70 + (c >> 16 & $ff)*$128e2f) >> 16
        If l0 & $8000
          ; midtone to highlight range
          l0 << 1
          l1 = ~l0
          Plot(x, y, ((l0*r2+l1*r1)>>24 & $ff)|((l0*g2+l1*g1)>>16 & $ff00)|((l0*b2+l1*b1)>>8 & $ff0000)|(c & $ff000000))
        Else
          ; shadow to midtone range
          l0 << 1
          l1 = ~l0
          Plot(x, y, ((l0*r1+l1*r0)>>24 & $ff)|((l0*g1+l1*g0)>>16 & $ff00)|((l0*b1+l1*b0)>>8 & $ff0000)|(c & $ff000000))
        EndIf
        x + 1
      Wend
      y + 1
    Wend
    StopDrawing()
  EndIf
  
EndProcedure

macros · Post by **macros** » Fri Feb 19, 2016 12:07 am

The speed is awesome!
Thank you for sharing it and thank you even more for explaining it with a plain code!

BasicallyPure · Post by **BasicallyPure** » Thu Feb 25, 2016 5:30 am

How about this procedure to produce quadtone images.
The 'shadow' parameter has been removed and is assumed to be zero.
C1 is the 'highlight' color and C2, C3, C4 are progressively darker mid-tones.
The transition points between colors, determined by luma, are at 50%, 25%, and 12.5%.
It seems to be the only options that work properly.
The default parameter values produce a high-contrast grayscale image.

Usually a color that is at maximum brightness (at least one RGB component = 255) works best for color parameters.

original grayscale images

after quadtone applied, C1,C2,C3,C4 = $00D9FF, $FFC981, $BC6BFF, $4948FF......................; C1,C2,C3,C4 = $60C9FF, $FFFF00, $FF00FF, $0000FF

Code: Select all

Procedure COLORIZE_QUADTONE(Image, C1=$FFFFFF, C2=$C0C0C0, C3=$808080, C4=$404040)
   
   Protected.i c, x, y, w, h
   Protected.i R1, G1, B1, R2, G2, B2, R3, G3, B3, R4, G4, B4
   Protected.i lum_1, lum_2
   
   ; expand color components
   R1 = Red(C1) : G1 = Green(C1) : B1 = Blue(C1)
   R2 = Red(C2) : G2 = Green(C2) : B2 = Blue(C2)
   R3 = Red(C3) : G3 = Green(C3) : B3 = Blue(C3)
   R4 = Red(C4) : G4 = Green(C4) : B4 = Blue(C4)
   
   ; process image
   If IsImage(Image) And StartDrawing(ImageOutput(Image))
         w = OutputWidth()
         h = OutputHeight()
         While y < h
            x = 0
            While x < w
               c = Point(x, y)
               
               ; calculate luma where lum_1 <= $1FFFF
               ; scale = R*0.299, G*0.587, B*0.114
               lum_1 = c & $FF * $266C0 : c >> 8 ; red
               lum_1 + c & $FF * $4B6E2 : c >> 8 ; grn
               lum_1 + c & $FF * $0EA63          ; blu
               lum_1 >> 10
               
               If lum_1 > $FFFF ; C1 to C2 range
                  lum_1 & $FFFF
                  lum_2 = ($FFFF - lum_1)>>1
                  Plot(x, y,
                       RGB((lum_1*R1 + lum_2*R2)>>16, ;red
                           (lum_1*G1 + lum_2*G2)>>16, ;grn
                           (lum_1*B1 + lum_2*B2)>>16));blu
                  
               ElseIf lum_1 > $7FFF ; C2 to C3 range
                  lum_1 & $7FFF
                  lum_2 = ($7FFF - lum_1)>>1
                  Plot(x, y,
                       RGB((lum_1*R2 + lum_2*R3)>>16, ;red
                           (lum_1*G2 + lum_2*G3)>>16, ;grn
                           (lum_1*B2 + lum_2*B3)>>16));blu
                  
               ElseIf lum_1 > $3FFF ; C3 to C4 range
                  lum_1 & $3FFF
                  lum_2 = ($3FFF - lum_1)>>1
                  Plot(x, y,
                       RGB((lum_1*R3 + lum_2*R4)>>16, ;red
                           (lum_1*G3 + lum_2*G4)>>16, ;grn
                           (lum_1*B3 + lum_2*B4)>>16));blu
               
               Else ; C4 to dark range
                  lum_1 & $1FFF
                  Plot(x, y,
                       RGB((lum_1*R4)>>16, ;red
                           (lum_1*G4)>>16, ;grn
                           (lum_1*B4)>>16));blu
               EndIf
               
               x + 1
            Wend
            y + 1
         Wend
         
      StopDrawing()
   EndIf
   
EndProcedure

Here is code that will calculate those pesky multipliers used in luma calculation.

Code: Select all

#maxLum  = $1FFFF ; value of lum when color = #White
#shift_n = 10     ; shift right this many bits at final calculation (>> n)

m.d = #maxLum / $FF * Pow(2,#shift_n)
m = Int(m)

; multipliers
Rm.i = m * 0.299
Gm.i = m * 0.587
Bm.i = m * 0.114

Debug "Rm = $" + RSet(Hex(Rm),6,"0")
Debug "Gm = $" + RSet(Hex(Gm),6,"0")
Debug "Bm = $" + RSet(Hex(Bm),6,"0")

Debug "Rm + Gm + Bm = $" + RSet(Hex(Rm+Gm+Bm),6,"0")

; verify the answers are correct
Debug "$" + Hex(((Rm+Gm+Bm)*255)) + " >> " + Str(#shift_n)
Debug "maxLum = $" + Hex((Rm*$FF + Gm*$FF + Bm*$FF)>>#shift_n) + " (want this to equal $" + Hex(#maxLum)

wilbert · Post by **wilbert** » Thu Feb 25, 2016 6:44 am

BasicallyPure wrote:How about this procedure to produce quadtone images.

The images you are showing look great !

I'm still working on it myself also.

I want to try 25%, 50% and 75% but maybe the points you chose are better.
While looking at this thread http://www.purebasic.fr/english/viewtop ... 17&t=64915 started by Keya and reading a bit more, I noticed the grayscale conversion functions used by professional software do the weighting formula in a non liner colorspace. The results look more natural but it takes more computing time so I'm experimenting with look up tables at the moment for grayscale conversion.

PureBasic Forums - English

Colorize module

Colorize module

Re: Colorize module

Re: Colorize module

Re: Colorize module

Re: Colorize module

Re: Colorize module