Intel® Fortran Compiler
Build applications that can scale for the future with optimized code designed for Intel® Xeon® and compatible processors.

Kmeans++

JohnNichols
Valued Contributor II
686 Views

I have been playing with KMEANS for a while, but I wanted to try it in Fortran.  

Anyway, here is a F90 version of the f77 version on the web.  Rosettacode.org.

 

 !**********************************************************************
    ! KMPP - K-Means++ - Traditional data clustering with a special initialization
    ! Public Domain - This program may be used by any person for any purpose.
    !
    ! Origin:
    !    Hugo Steinhaus, 1956
    !
    ! Refer to:
    !    "kmeans++: the advantages of careful seeding"
    !    David Arthur and Sergei Vassilvitskii
    !    Proceedings of the eighteenth annual ACM-SIAM symposium
    !      on Discrete algorithms, 2007
    !
    !____Variable_______I/O_______Description___________________Type_______
    !    X(P,N)         In        Data points                   Real
    !    P              In        Dimension of the data         Integer
    !    N              In        Number of points              Integer
    !    K              In        # clusters                    Integer
    !    C(P,K)         Out       Center points of clusters     Real
    !    Z(N)           Out       What cluster a point is in    Integer
    !    WORK(N)        Neither                                 Real
    !    IFAULT         Out       Error code                    Integer
    !***********************************************************************
    SUBROUTINE KMPP (X, P, N, K, C, Z, WORK, IFAULT)

    IMPLICIT NONE
    INTEGER P, N, K, Z, IFAULT
    REAL X, C, WORK
    DIMENSION X(P,N), C(P,K), Z(N), WORK(N)

    !               constants
    INTEGER ITER                 ! maximum iterations
    REAL BIG                     ! arbitrary large number
    PARAMETER (ITER = 1000, BIG = 1E33)

    !                local variables
    INTEGER    H          ! count iterations
    INTEGER         I          ! count points
    INTEGER          I1         ! point marked as initial center
    INTEGER         J         ! count dimensions
    INTEGER         L          ! count clusters
    INTEGER         L0         ! present cluster ID
    INTEGER         L1          ! new cluster ID

    REAL  BEST                 ! shortest distance to a center
    REAL      D2                   ! squared distance
    REAL      TOT                  ! a total
    REAL      W                     ! temp scalar

    LOGICAL CHANGE             ! whether any points have been reassigned

    !***********************************************************************
    !           Begin.
    !***********************************************************************
    IFAULT = 0
    IF (K < 1 .OR. K > N) THEN       ! K out of bounds
        IFAULT = 3
        RETURN
    END IF
    DO I = 1, N                       ! clear Z
        Z(I) = 0
    END DO

    !***********************************************************************
    !        initial centers
    !***********************************************************************
    DO I = 1, N
        WORK(I) = BIG
    END DO

    CALL RANDOM_NUMBER (W)
    I1 = MIN(INT(W * FLOAT(N)) + 1, N)  ! choose first center at random
    DO J = 1, P
        C(J,1) = X(J,I1)
    END DO

    DO L = 2, K                    ! initialize other centers
        TOT = 0.
        DO I = 1, N                     ! measure from each point
            BEST = WORK(I)
            D2 = 0.                         ! to prior center
            DO J = 1, P
                D2 = D2 + (X(J,I) - C(J,L-1)) **2  ! Squared Euclidean distance
                IF (D2 .GE. BEST) GO TO 10               ! needless to add to D2
            END DO                          ! next J
            IF (D2 < BEST) BEST = D2          ! shortest squared distance
            WORK(I) = BEST
10          TOT = TOT + BEST             ! cumulative squared distance
        END DO                      ! next data point

        !***********************************************************************
        ! Choose center with probability proportional to its squared distance
        !     from existing centers.
        !***********************************************************************
        CALL RANDOM_NUMBER (W)
        W = W * TOT    ! uniform at random over cumulative distance
        TOT = 0.
        DO I = 1, N
            I1 = I
            TOT = TOT + WORK(I)
            IF (TOT > W) GO TO 20
        END DO                ! next I
20      CONTINUE
        DO J = 1, P         ! assign center
            C(J,L) = X(J,I1)
        END DO
    END DO               ! next center to initialize

    !***********************************************************************
    !                      main loop
    !***********************************************************************
    DO H = 1, ITER
        CHANGE = .FALSE.

        !             find nearest center for each point
        DO I = 1, N
            L0 = Z(I)
            L1 = 0
            BEST = BIG
            DO L = 1, K
                D2 = 0.
                DO J = 1, P
                    D2 = D2 + (X(J,I) - C(J,L)) **2
                    IF (D2 .GE. BEST) GO TO 30
                END DO
30              CONTINUE
                IF (D2 < BEST) THEN           ! new nearest center
                    BEST = D2
                    L1 = L
                END IF
            END DO        ! next L

            IF (L0 .NE. L1) THEN
                Z(I) = L1                   !  reassign point
                CHANGE = .TRUE.
            END IF
        END DO         ! next I
        IF (.NOT. CHANGE) RETURN      ! success

        !***********************************************************************
        !           find cluster centers
        !***********************************************************************
        DO L = 1, K              ! zero population
            WORK(L) = 0.
        END DO
        DO L = 1, K               ! zero centers
            DO J = 1, P
                C(J,L) = 0.
            END DO
        END DO

        DO I = 1, N
            L = Z(I)
            WORK(L) = WORK(L) + 1.             ! count
            DO J = 1, P
                C(J,L) = C(J,L) + X(J,I)         ! add
            END DO
        END DO

        DO L = 1, K
            IF (WORK(L) < 0.5) THEN          ! empty cluster check
                IFAULT = 1                     ! fatal error
                RETURN
            END IF
            W = 1. / WORK(L)
            DO J = 1, P
                C(J,L) = C(J,L) * W     ! multiplication is faster than division
            END DO
        END DO

    END DO                   ! next H
    IFAULT = 2                ! too many iterations
    RETURN

    END  ! of KMPP


    !***********************************************************************
    !             test program (extra credit #1)
    !***********************************************************************
    PROGRAM TPEC1
    IMPLICIT NONE
    Integer, PARAMETER :: N = 30, P = 2, K = 6
    REal, Parameter :: TWOPI = 6.2831853
    INTEGER I, L, Z(N), IFAULT
    REAL X(P,N), C(P,K), R, THETA, W, WORK(N)

    !             Begin
    CALL RANDOM_SEED()
    DO I = 1, N                      ! random points over unit circle
        CALL RANDOM_NUMBER (W)
        R = SQRT(W)                      ! radius
        CALL RANDOM_NUMBER (W)
        THETA = W * TWOPI                ! angle
        write(*,*)"A   :: ",I, R,THETA
        X(1,I) = R * COS(THETA)          ! Cartesian coordinates
        X(2,I) = R * SIN(THETA)
        !write(*,*)I,X(1,I),X(2,I)
    END DO

    !             Call subroutine
    CALL KMPP (X, P, N, K, C, Z, WORK, IFAULT)
    PRINT *, 'kmpp returns with error code ', IFAULT

    !            Print lists of points in each cluster
    DO L = 1, K
        PRINT *, 'Cluster ', L, ' contains points: '
10      FORMAT (I6, $)
20      FORMAT ()
        DO I = 1, N
            IF (Z(I) .EQ. L) PRINT 10, I
        END DO
        PRINT 20
    END DO

    !         Write CSV file with Y-coordinates in different columns by cluster
    OPEN (UNIT=1, FILE='tpec1.csv', STATUS='UNKNOWN', IOSTAT=IFAULT)
    IF (IFAULT .NE. 0) PRINT *, 'tpec1: trouble opening file'
30  FORMAT (F8.4, $)
40  FORMAT (',', $)
50  FORMAT (F8.4)
    DO I = 1, N
        WRITE (UNIT=1, FMT=30, IOSTAT=IFAULT) X(1,I)
        IF (IFAULT .NE. 0) PRINT *, 'tpec1: trouble writing X-coord'
        DO L = 1, Z(I)                     ! one comma per cluster ID
            WRITE (UNIT=1, FMT=40, IOSTAT=IFAULT)
            IF (IFAULT .NE. 0) PRINT *, 'tpec1: trouble writing comma'
        END DO
        WRITE (UNIT=1, FMT=50, IOSTAT=IFAULT) X(2,I)
        IF (IFAULT .NE. 0) PRINT *, 'tpec1: trouble writing Y-coord'
    END DO

    !           Write the centroids in the far column
    DO L = 1, K
        WRITE (UNIT=1, FMT=30, IOSTAT=IFAULT) C(1,L)
        IF (IFAULT .NE. 0) PRINT *, 'tpec1: trouble writing X-coord'
        DO I = 1, K+1
            WRITE (UNIT=1, FMT=40, IOSTAT=IFAULT)
            IF (IFAULT .NE. 0) PRINT *, 'tpec1: trouble writing comma'
        END DO
        WRITE (UNIT=1, FMT=50, IOSTAT=IFAULT) C(2,L)
        IF (IFAULT .NE. 0) PRINT *, 'tpec1: trouble writing Y-coord'
    END DO
    CLOSE (UNIT=1)

    END  ! of test program

 

 

it needs more work, but it is a neat little program.  

0 Kudos
8 Replies
JohnNichols
Valued Contributor II
634 Views
n = RGBnode%nn
    call RGBnode%fill(1,n )

I have to declare an integer n(2) and set it to pass n RGBnode%nn into the subroutine, it will not pass directly - I just do not understand why I need to do this?

 

Line 194 is the relevant subroutine. 

 module rgbimage_m

    implicit none

    private
    public :: rgbimage
    INTEGER, PARAMETER :: dp = selected_real_kind(15, 307)
    REAL (KIND=dp) :: g = 9.806, pi = 3.14159265D0


    type rgbimage
        !! usage
        !!    1) init
        !!    2a) fill_image
        !!      or
        !!    2b) set_pixel
        !!    3) normalize
        !!    4) write

        integer Colred(3), colgre(3), colblu(3), Colyel(3)
        
        integer, dimension(:,:,:), allocatable :: input

        integer, dimension(:,:,:), allocatable :: rgb
        !! pixel arrays of rgb values
        !! indices (i,j,k)
        !!    i: position x_i
        !!    j: position y_j
        !!    k=1: red, k=2: green, k=3: blue

        integer :: n(2) 
        !! image dimensions: [height, width]
        integer :: nn(2) 
        !! image input dimensions: [height, width]



    contains
    procedure :: init       => rgbimage_init              ! inits image
    procedure :: fill_image => rgbimage_fill_image        ! fill image with constant rgb value
    procedure :: get_pixel  => rgbimage_get_pixel         ! gets one pixel
    procedure :: normalize  => rgbimage_normalize         ! normalizes all pixels onto range [0, 255]
    procedure :: set_pixel  => rgbimage_set_pixel         ! sets one pixel
    procedure :: write      => rgbimage_write             ! outputs image to file
    procedure :: fill      => rgbimage_fill               ! outputs image to file

    procedure, private :: inside => rgbimage_inside
    procedure, private :: valid  => rgbimage_valid
    end type

    contains



    subroutine rgbimage_init(this, height, width)
    !! initialize image.
    !! sets dimensions, allocates pixels and sets colors to 0.

    class(rgbimage), intent(out) :: this
    integer,         intent(in)  :: height, width
    integer astat
    
     integer, parameter :: ihpixf =  4096, jvpixf = 4096 ! pixel size, jvpixf

    this%n = [height, width]
    
    this%nn = [ihpixf,jvpixf]
    allocate (this%rgb(height,width,3), source=0, STAT = astat)
    if(astat == 0) then
        write(*,*)"Allocated matrix"
    else
        stop "Failed to allocate matrix"
    end if
    allocate(this%input(3,ihpixf,jvpixf),STAT = astat)
    
    if(astat == 0) then
        write(*,*)"Allocated matrix"
    else
        stop "Failed to allocate matrix"
    end if
    this%Colgre(1) = 248
    this%Colgre(2) = 248
    this%Colgre(3) = 255

    this%Colyel(1) = 248
    this%Colyel(2) = 100
    this%Colyel(3) = 100

    this%Colblu(1) = 50
    this%Colblu(2) = 56
    this%Colblu(3) = 200

    this%Colred(1) = 200
    this%Colred(2) = 256
    this%Colred(3) = 256
    return
    end subroutine

    logical function rgbimage_valid(this, check_rgb_vals)
    !! checks if the image has valid dimensions and optionally valid rgb values.

    class(rgbimage), intent(in)           :: this
    logical,         intent(in), optional :: check_rgb_vals
    !! check if rgb values are in allowed range [0, 255]?
    !! default: dont check

    ! always check that dimensions match
    rgbimage_valid = ( all(this%n > 0)                     .and. &
        &               (size(this%rgb, dim=1) == this%n(1)) .and. &
        &               (size(this%rgb, dim=2) == this%n(2)) .and. &
        &               (size(this%rgb, dim=3) == 3)               )

    ! optionally: check if rgb values are in allowed range
    if (present(check_rgb_vals)) then
        if (check_rgb_vals) rgbimage_valid = ( rgbimage_valid       .and. &
            &                                   (all(this%rgb >= 0))  .and. &
            &                                   (all(this%rgb <= 255))      )
    end if


    end function

    logical function rgbimage_inside(this, x, y)
    !! checks if given coordinates are inside the image

    class(rgbimage), intent(in) :: this
    integer,         intent(in) :: x, y

    rgbimage_inside = ((x > 0) .and. (x <= this%n(1)) .and. (y > 0) .and. (y <= this%n(2)))
    end function

    subroutine rgbimage_set_pixel(this, x, y, rgb)
    class(rgbimage), intent(inout) :: this
    integer,         intent(in)    :: x, y
    !! coordinates
    integer,         intent(in)    :: rgb(3)
    !! red, green, blue values

    if (this%inside(x, y)) then
        ! use given data at first
        this%rgb(x,y,:) = rgb

        ! check if given data was out of bounds
        where     (this%rgb(x,y,:) > 255)
            this%rgb(x,y,:) = 255
        elsewhere (this%rgb(x,y,:) < 0)
            this%rgb(x,y,:) = 0
        end where
    end if
    end subroutine

    function rgbimage_get_pixel(this, x, y) result(rgb)
    class(rgbimage), intent(in) :: this
    integer,         intent(in) :: x, y
    !! coordinates
    integer                     :: rgb(3)
    !! red, green, blue values

    if (this%inside(x, y)) then
        rgb = this%rgb(x,y,:)
    else
        rgb = 0
    end if
    end function

    subroutine rgbimage_normalize(this)
    !! normalize colors to be in range [0, 255]

    class(rgbimage), intent(inout) :: this

    where     (this%rgb(:,:,:) > 255)
        this%rgb(:,:,:) = 255
    elsewhere (this%rgb(:,:,:) < 0)
        this%rgb(:,:,:) = 0
    end where
    end subroutine

    subroutine rgbimage_fill_image(this, rgb)
    !! fill whole image with given rgb values.

    class(rgbimage), intent(inout) :: this
    integer,         intent(in)    :: rgb(3)
    !! red, green, blue values

    integer :: i

    if (this%valid()) then
        do i = 1, 3
            this%rgb(:,:,i) = rgb(i)
        end do
    end if
    end subroutine
    
    subroutine rgbimage_fill(this, i, n)
    implicit none
    
    class(rgbimage), intent(out) :: this
    integer, intent(in) :: i, n(2)
    
    integer j
    
    do 100 j = 1, n(1)
    write(*,*)"Here  :: ",i,j
100 end do    
    
    
    
    
    return 
    end subroutine
    
    

    subroutine rgbimage_write(this, fname)
    class(rgbimage), intent(in) :: this
    character(*),    intent(in) :: fname
    character*100 :: error
    !! file path, e.g. "tmp/out.ppm"

    integer :: iounit, ios, i,j,k

    open (newunit=iounit, file=fname, iostat=ios, action='WRITE')
    if (ios /= 0) then
        write(error,100)fname
100     format(    "Error opening file: ", A20)
        stop 'Error opening file'
    endif

    ! write header
    write (iounit, '(A)')         'P6'
    write (iounit, '(I0, A, I0)') this%n(1), " ", this%n(2)
    write (iounit, '(A)')         '255'

    do i = 1, this%n(1)
        do j = 1, this%n(2)
            write (iounit, '(3A1)', advance='no') [(achar(this%rgb(i,j,k)), k=1,3)]
        end do
    end do

    close (unit=iounit, iostat=ios)
    if (ios /= 0) error stop "Error closing file"
    end subroutine


    ! --------------------------------------------
    !
    ! Notes
    ! o With a parameter ipixout set at 1, 2 or others,
    !   this subroutine will generate PPM-P6(binary), PPM-P3(text),
    !   or BMP(24bit depth without color table).
    !
    ! o Some parts follow DEC-FORTRAN convention that had been defacto-standard long ago.
    !   Some compilers today may not accept if "ipixout" is set other than 2.
    !
    ! o g77 (ver. 3.3.3) works for all three choices.
    ! o Intel compiler (ver. 9 or later) works for all three choices.
    !
    ! --------------------------------------------
    !
    subroutine pixout(rgb,nframe)
    implicit none
    ! interface arg.
    integer ihpixf, jvpixf
    parameter(ihpixf = 4096, jvpixf = 4096) ! pixel size, eacg must be multiple of 4, if BMP is chosen as output format.
    character*1 rgb(3,ihpixf,jvpixf)      ! RGB data array
    integer nframe
    ! local
    character*12 fnameout
    integer i, j, k
    integer itmp, icnt
    character*14 frmtstr
    character*54 headmsw
    character*4  byt4
    character*2  byt2
    ! choices
    integer ipixout
    parameter(ipixout = 3) ! 1 / 2 / other= PPM6, PPM3, BMP(24bit)

    if (ipixout .EQ. 1) then

        ! PPM P6

        write(fnameout,'(''smpl'',i3.3,''.ppm'')') nframe ! name of PPM file
        open(unit=2,file=fnameout,status='unknown')
        write(*,*) 'Now writing PPM (P6) file : ', fnameout
        ! header
        write(2,'(''P6'', 2(1x,i4),'' 255 '',$)')ihpixf, jvpixf
        ! image data
        itmp = ihpixf ! jvpixf ! 3
        write(frmtstr,'(''('',i8.8,''A,$)'')') itmp     ! make output "format"
        write(2,fmt=frmtstr) (((rgb(k,i,j),k=1,3),i=1,ihpixf),j=jvpixf,1,-1)

        close(2)

    else if (ipixout .EQ. 2) then

        ! PPM P3 ! rather "safer" choice for many Fortran compiler(s).

        write(fnameout,'(''smpl'',i3.3,''.ppm'')') nframe ! name of PPM file
        open(unit=2,file=fnameout,status='unknown')
        write(*,*) 'Now writing PPM (P3) file : ', fnameout
        ! header
        write(2,'(A)') 'P3'
        write(2,'(2(1x,i4),'' 255 '')')  ihpixf, jvpixf
        icnt = 0
        ! image data
        do j = jvpixf, 1, -1                              ! here, j (vertical address) runs from top to bottom.
            do i = 1, ihpixf, 1
                do k = 1, 3
                    itmp = ichar(rgb(k,i,j))
                    icnt = icnt + 4
                    if (icnt .LT. 60) then
                        write(2,fmt='(1x,i3,$)') itmp                 ! "$" is not standard.
                    else
                        write(2,fmt='(1x,i3)') itmp
                        icnt = 0
                    endif
                enddo
            enddo
        enddo
        write(2,'(A)') ' '
        close(2)

    else

        ! BMP (24bit depth)... this part works only when width is multiple of 4.

        itmp = mod(ihpixf, 4)
        if (itmp .NE. 0) then
            write(*,*) 'width must be multiple of 4'
            stop
        endif

        write(fnameout,'(''smpl'',i3.3,''.bmp'')') nframe ! name of BMP file
        open(unit=2,file=fnameout,status='unknown')
        write(*,*) 'Now writing BMP(24bit) file : ', fnameout
        ! header 1 (file header ; 1--14 byte)
        headmsw( 1: 2) = 'BM'             ! declaring this is BMP file
        itmp = 54 + ihpixf * jvpixf * 3 ! total file size = header + data
        call num2bit4(itmp,byt4)
        headmsw( 3: 6) = byt4(1:4)
        itmp = 0                        ! may be 0
        call num2bit2(itmp,byt2)
        headmsw( 7:  = byt2(1:2)
        itmp = 0                        ! may be 0
        call num2bit2(itmp,byt2)
        headmsw( 9:10) = byt2(1:2)
        itmp = 54                       ! must be 54 : total length of header
        call num2bit4(itmp,byt4)
        headmsw(11:14) = byt4(1:4)
        ! header 2 (bit-map header ; 13--54 byte)
        itmp = 40                       ! must be 40 : length of bit-map header
        call num2bit4(itmp,byt4)
        headmsw(15:18) = byt4(1:4)
        itmp = ihpixf                   ! width
        call num2bit4(itmp,byt4)
        headmsw(19:22) = byt4(1:4)
        itmp = jvpixf                   ! height
        call num2bit4(itmp,byt4)
        headmsw(23:26) = byt4(1:4)
        itmp = 1                        ! must be 1
        call num2bit2(itmp,byt2)
        headmsw(27:28) = byt2(1:2)
        itmp = 24                       ! must be 24 : color depth in bit.
        call num2bit2(itmp,byt2)
        headmsw(29:30) = byt2(1:2)
        itmp = 0                        ! may be 0 : compression method index
        call num2bit4(itmp,byt4)
        headmsw(31:34) = byt4(1:4)
        itmp = 0                        ! may be 0 : file size if compressed
        call num2bit4(itmp,byt4)
        headmsw(35:38) = byt4(1:4)
        itmp = 0                        ! arbit. : pixel per meter, horizontal
        call num2bit4(itmp,byt4)
        headmsw(39:42) = byt4(1:4)
        itmp = 0                        ! arbit. : pixel per meter, vertical
        call num2bit4(itmp,byt4)
        headmsw(43:46) = byt4(1:4)
        itmp = 0                        ! may be 0 here : num. of color used
        call num2bit4(itmp,byt4)
        headmsw(47:50) = byt4(1:4)
        itmp = 0                        ! may be 0 here : num. of important color
        call num2bit4(itmp,byt4)
        headmsw(51:54) = byt4(1:4)

        ! writing header part
        write(2,'(a54,$)') headmsw(1:54)
        ! image data
        itmp = ihpixf * jvpixf * 3
        write(frmtstr,'(''('',i8.8,''A,$)'')') itmp
        write(2,fmt=frmtstr) (((rgb(k,i,j),k=3,1,-1),i=1,ihpixf),j=1,jvpixf) ! writing in BGR order, not RGB.
        close(2)

    endif

    return
    end subroutine pixout

    ! --------------------------------------
    ! convert integer values to 4 8-bit characters
    ! --------------------------------------

    subroutine num2bit4(inum,byt4)
    implicit none
    integer inum
    character*4 byt4
    integer itmp1, itmp2
    itmp1 = inum
    itmp2 = itmp1 / 256**3
    byt4(4:4) = char(itmp2)
    itmp1 =-itmp2 * 256**3 +itmp1
    itmp2 = itmp1 / 256**2
    byt4(3:3) = char(itmp2)
    itmp1 =-itmp2 * 256**2 +itmp1
    itmp2 = itmp1 / 256
    byt4(2:2) = char(itmp2)
    itmp1 =-itmp2 * 256    +itmp1
    byt4(1:1) = char(itmp1)
    return
    end subroutine num2bit4

    ! --------------------------------------
    ! convert integer values to 2 8-bit characters
    ! --------------------------------------

    subroutine num2bit2(inum,byt2)
    implicit none
    integer inum
    character*2 byt2
    integer itmp1, itmp2
    itmp1 = inum
    itmp2 = itmp1 / 256
    byt2(2:2) = char(itmp2)
    itmp1 =-itmp2 * 256 + itmp1
    byt2(1:1) = char(itmp1)
    return
    end subroutine num2bit2


    ! --------------------------------------------
    !   fill rgb data array with something
    ! --------------------------------------------


    ! --------------------------------------------
    ! Sample FORTRAN 77 program to make PPM / BMP
    !
    !  1998-Apr-03 (long long ago)
    !      :
    !  2005-Jul-22 (added P3)
    !      :
    !  2006 Dec-22 (added BMP)
    !
    !  Image array rgb(3,!,!) is filled in subroutine mkbitmap()
    !  and
    !  saved in PPM or BMP format in subroutine pixout().
    !
    !                                   K. Hayashi
    ! --------------------------------------------
    !
    subroutine pixelout(input, mag)
    implicit none
    integer ihpixf, jvpixf
    parameter(ihpixf = 4096, jvpixf = 4096) ! pixel size
    character*1 rgb(3,ihpixf,jvpixf) ! RGB image array
    integer nframe, nf2, i, j
    integer input(3, ihpixf,jvpixf)
    integer Colred(3)
    real mag



    do nframe = 1, 1
        nf2 = int (mag)
        call mkbitmap(rgb,nf2, input)
        call pixout(rgb,nf2)
    enddo


    end subroutine pixelout



    subroutine mkbitmap(rgb,nframe, input)
    implicit none
    integer nframe
    integer ihpixf, jvpixf
    parameter(ihpixf = 4096, jvpixf = 4096) ! pixel size
    character*1 rgb(3,ihpixf,jvpixf) !      RGB pixel data array
    integer input(3, ihpixf,jvpixf)
    ! local
    real*8  red, gre, blu
    integer ired, igre, iblu
    real*8  ofst
    parameter(ofst = 0.7D+00)
    integer i, j, itmp
    real*8   aa, bb, cc, rr, xx, yy, tt
    integer ichoice
    parameter(ichoice = 0) ! .... choice
    real*8  pi
    parameter(pi = 3.14159265358979D+00)


    if (ichoice .EQ. 0) then
        do 100 j = 1, jvpixf
            do 100 i = 1, ihpixf

                rgb(1,i,j) = char(input(1,i,j)) ! red
                rgb(2,i,j) = char(input(2,i,j)) ! green
                rgb(3,i,j) = char(input(3,i,j)) ! blue
100     continue

    endif

    return
    end subroutine mkbitmap


    ! --------------------------------------------
    ! end of this file, thank you.
    ! --------------------------------------------
    end module

   

jimdempseyatthecove
Black Belt
628 Views

I haven't tried your code. This seems like a regression error (bug) where the UDT%var as argument is seen as UDT alone.

Try: call RGBnode%fill(1,(RGBnode%nn))

The extra parenthesis will construct a copy of the member variable (unfortunate as it is, but removes the necessity for a user defined temp.

Jim Dempsey

JohnNichols
Valued Contributor II
610 Views
jimdempseyatthecove
Black Belt
576 Views

John, et. al.

I must caution that that fix only applies to arguments that are intended to be changed/defined by the called procedure.

Jim Dempsey

JohnNichols
Valued Contributor II
534 Views

 

In pseudo-code
Module AA
TYPE A 
INTEGER N(2)
End TYPE
FUNCTION START(N)
set N(2)
end function
FUNCTION USE(N)
USE N for work
end function
END Module AA

MAIN PROGRAM
TYPE A == set new
CALL A%Start(A%N)
!! A%N has now two integer values
CALL %Use((A%N))
!! Do not change it -- it was set in Start.

End main

 

This is how it is working at the moment, is that the intended use of the (())

 

JohnNichols
Valued Contributor II
510 Views

KMeans++ ResultKMeans++ Result

It has been fun to get the KMeans working in Intel Fortran and putting a drawing program with it, so I get straight plots.  

 

cean
New Contributor I
433 Views

Nice. Where I can get more information about these colors in rgb.f90? I want to make a color scale, from blue to red, 6~8 values in between. Thanks.

 

        case(1)                         !Medium Violet Red
            this%colors(i,1) = 199
            this%colors(i,2) = 21
            this%colors(i,3) = 133
            this%list(i) = "MVP"

 

JohnNichols
Valued Contributor II
407 Views

In the ZIP file are three solutions in Fortran.  They work on VS 2019 - but PSXE 2020.4.  I had to go back to get a working Fortran with the latest VS.  

The KMEANS solution works on a simple set of FFT data.  The data is enclosed.  

RGB.f90 is enclosed in the KM++ solution folder.  It has the 140 colors of the WEB Safe W3 coded into the routine.  You call a color by number.   ignore the RGB in the Kmeans it is an old version. 

I enclose the web safe PDF that lists all the colors.  I have the subroutine set up for 140 colors, but only 70 have the full RGB codes in place. If you see Pink it is a place holder.   I have found it is much easier to stick to the color names in this list.  

I include the code so you have samples of the calling routines.  

This is not production code, it is for experimenting with data. 

Reply