Path: news.mathworks.com!not-for-mail
From: <HIDDEN>
Newsgroups: comp.soft-sys.matlab
Subject: Re: mxGetNumberOfDimensions missing ?
Date: Thu, 5 Nov 2009 09:31:01 +0000 (UTC)
Organization: Boeing
Lines: 359
Message-ID: <hcu60l$n6v$1@fred.mathworks.com>
References: <2009102112385016807-spammersgohere@spaminvalid> <hbnv28$c5h$1@fred.mathworks.com> <2009102115045916807-spammersgohere@spaminvalid> <hbo6a4$k5l$1@fred.mathworks.com> <2009102311375416807-spammersgohere@spaminvalid> <hc56i5$3p4$1@fred.mathworks.com> <hc5jrq$b1n$1@fred.mathworks.com> <hc7ebk$8nl$1@fred.mathworks.com> <2009110417351816807-spammersgohere@spaminvalid> <hctanf$dbk$1@fred.mathworks.com> <2009110421001016807-spammersgohere@spaminvalid>
Reply-To: <HIDDEN>
NNTP-Posting-Host: webapp-03-blr.mathworks.com
Content-Type: text/plain; charset="ISO-8859-1"
Content-Transfer-Encoding: 8bit
X-Trace: fred.mathworks.com 1257413461 23775 172.30.248.38 (5 Nov 2009 09:31:01 GMT)
X-Complaints-To: news@mathworks.com
NNTP-Posting-Date: Thu, 5 Nov 2009 09:31:01 +0000 (UTC)
X-Newsreader: MATLAB Central Newsreader 756104
Xref: news.mathworks.com comp.soft-sys.matlab:582661


Geico Caveman <spammers-go-here@spam.invalid> wrote in message <2009110421001016807-spammersgohere@spaminvalid>...
> On 2009-11-04 18:45:19 -0700, "James Tursa" 
> <aclassyguy_with_a_k_not_a_c@hotmail.com> said:
> 
> Can I email it to you ?
> 

Yes. Just replace the c with a k in the first word of the e-mail address above.

> Sounds interesting. Maybe you could post this to Matlab's code exchange 
> once it is ready.

Will do. I will let you know when it gets posted to the FEX. In the meantime, see below ...

> Yes. Basically  have three indices. I would want to be able to access 
> the 3D array pointed to for manipulation in Matlab and Fortran.

OK, here is a dramatically pared down example of my code. I clipped out everything except the 3D stuff. It is a simple mex function that performs the calculation sin(A)+7 where A must be a 3D double input. Uses Fortran pointers, accesses the whole arrays through the pointers, accesses the dimensions directly, and allocates memory for a copy and makes a copy and then deallocates it. Enjoy.

James Tursa

!----------------------------------------------------------------------------------------------------------

!-----------------------------------------------------------------------------------
! Programmer: James Tursa
! Function to demonstrate use of Fortran pointers on a 3D double input.
! Performs the calculation B = sin(A) + 7
! Intended for very small dimension sizes because of the printing
!

#include "fintrf.h"

!\
! The following macros are needed for older versions of MATLAB that do not have these
! macros defined in the fintrf.h file.
!/

#ifndef mwSize
#define mwSize integer(4)
#endif

#ifndef mwPointer
#define mwPointer integer(4)
#endif

#ifndef mwIndex
#define mwIndex integer(4)
#endif

!---------------------------------------------------------------------
!---------------------------------------------------------------------
!---------------------------------------------------------------------

      module MatlabAPImex
      
!-----------------------------------------------------      
! Interface definitions for MATLAB API mex functions
!-----------------------------------------------------
      
      interface
!-----
      integer(4) function mexPrintf(message)
      character(len=*), intent(in) :: message
      end function mexPrintf
!-----
      end interface

      end module MatlabAPImex

!---------------------------------------------------------------------
!---------------------------------------------------------------------
!---------------------------------------------------------------------

      module MatlabAPImx
      
!-----------------------------------------------------      
! Interface definitions for MATLAB API mx functions
!-----------------------------------------------------
      
      interface
!-----
      mwPointer function mxDuplicateArray(pm)
      mwPointer, intent(in) :: pm
      end function mxDuplicateArray
!-----
      subroutine mxFree(ptr)
      mwPointer, intent(in) :: ptr
      end subroutine mxFree
!-----
      mwPointer function mxGetDimensions(pm)
      mwPointer, intent(in) :: pm
      end function mxGetDimensions
!-----
      mwSize function mxGetNumberOfDimensions( mx )
      mwPointer, intent(in) :: mx
      end function mxGetNumberOfDimensions
!-----
      mwPointer function mxGetPr( mx )
      mwPointer, intent(in) :: mx
      end function mxGetPr
!-----
      integer(4) function mxIsDouble( mx )
      mwPointer, intent(in) :: mx
      end function mxIsDouble
!-----
      integer(4) function mxIsSparse(pm)
      mwPointer, intent(in) :: pm
      end function mxIsSparse
!-----
      mwPointer function mxMalloc(n)
      mwSize, intent(in) :: n
      end function mxMalloc
!-----
      end interface

      end module MatlabAPImx

!---------------------------------------------------------------------
!---------------------------------------------------------------------
!---------------------------------------------------------------------

      module MatlabAPIfp
      use MatlabAPImx
      
!-----------------------------------------------------      
! Interface definitions for Fortan Pointer functions
!-----------------------------------------------------
      
      interface fpGetPr3
          module procedure fpGetPr3Double
      end interface
      
      interface fpAllocate
          module procedure fpAllocate3Double
      end interface
      
      interface fpDeallocate
          module procedure fpDeallocate3Double
      end interface
      
!-----------------------------------------------------      

      contains
      
!-----------------------------------------------------      
! Specific Fortran Pointer functions
!-----------------------------------------------------
      
!----------------------------------------------------------------------
      function fpGetPr3Double( mx ) result(fp)
      implicit none
      real(8), pointer :: fp(:,:,:)
!-ARG
      mwPointer, intent(in) :: mx
!-COM
      real(8), pointer :: Apx3(:,:,:)
      common /MatlabAPI_COMA/ Apx3
!-LOC
      mwPointer :: pr
      mwSize :: ndim
      mwSize, pointer :: dims(:)
      mwSize :: dimz(3)
!-----
      ndim = mxGetNumberOfDimensions(mx)
      if( mxIsDouble(mx) == 1 .and. mxIsSparse(mx) == 0 .and.           &
     &    ndim <= 3 ) then
          pr = mxGetPr( mx )
          dims => fpGetDimensions( mx )
          dimz = 1
          dimz(1:ndim) = dims
          call MatlabAPI_COM_Apx3( %val(pr), dimz )
          fp => Apx3
      else
          nullify( fp )
      endif
      return
      end function fpGetPr3Double
!----------------------------------------------------------------------
      function fpGetDimensions( mx ) result(fp)
      implicit none
      mwSize, pointer :: fp(:)
!-ARG
      mwPointer, intent(in) :: mx
!-COM
      mwSize, pointer :: Dpx(:)
      common /MatlabAPI_COMD/ Dpx
!-LOC
      mwSize :: ndim
      mwPointer :: dims
!-----
      ndim = mxGetNumberOfDimensions( mx )
      dims = mxGetDimensions( mx )
      call MatlabAPI_COM_Dpx(%val(dims), ndim)
      fp => Dpx
      end function fpGetDimensions
!----------------------------------------------------------------------
      function fpAllocate3Double( n1, n2, n3 ) result(fp)
      implicit none
      real(8), pointer :: fp(:,:,:)
!-ARG
      mwSize, intent(in) :: n1, n2, n3
!-COM
      real(8), pointer :: Apx3(:,:,:)
      common /MatlabAPI_COMA/ Apx3
!-LOC
      mwPointer :: mxmemory
!-----
      mxmemory = mxMalloc(n1*n2*n3*8)
      call MatlabAPI_COM_Apx3( %val(mxmemory), (/n1,n2,n3/) )
      fp => Apx3
      return
      end function fpAllocate3Double
!----------------------------------------------------------------------
      subroutine fpDeallocate3Double( fp )
      implicit none
!-ARG
      real(8), pointer, intent(inout) :: fp(:,:,:)
!-LOC
      mwPointer :: mxmemory
!-----
      if( associated(fp) ) then
          mxmemory = %loc(fp)
          call mxFree(mxmemory)
          nullify(fp)
      endif
      return
      end subroutine fpDeallocate3Double

      end module MatlabAPIfp

!---------------------------------------------------------------------
!---------------------------------------------------------------------
!---------------------------------------------------------------------
      
!----------------------------------------------------------------------
! Specific Fortan Helper functions. Not contained in the module because
! we need an implicit interface to get the %val() construct to work
! properly in the calling routine. Passing the appropriate pointer back
! in a COMMON block. Looks awkward, but works beautifully.
!----------------------------------------------------------------------
      subroutine MatlabAPI_COM_Apx3( A, DIMS )
      implicit none
!-ARG
      mwSize, intent(in) :: DIMS(3)
      real(8), target, intent(in) :: A(DIMS(1),DIMS(2),DIMS(3))
!-COM
      real(8), pointer :: Apx3(:,:,:)
      common /MatlabAPI_COMA/ Apx3
!-----
      Apx3 => A
      return
      end subroutine MatlabAPI_COM_Apx3
!----------------------------------------------------------------------      
      subroutine MatlabAPI_COM_Dpx(dims, ndim)
      implicit none
!-ARG
      mwSize, intent(in) :: ndim
      mwSize, target, intent(in) :: dims(ndim)
!-COM
      mwSize, pointer :: Dpx(:)
      common /MatlabAPI_COMD/ Dpx
!-----
      Dpx => dims
      return
      end subroutine MatlabAPI_COM_Dpx

!---------------------------------------------------------------------
!---------------------------------------------------------------------
!---------------------------------------------------------------------

      subroutine mexFunction(nlhs, plhs, nrhs, prhs)
      use MatlabAPImex
      use MatlabAPImx
      use MatlabAPIfp
      implicit none
!-ARG
      mwPointer plhs(*), prhs(*)
      integer*4 nlhs, nrhs
!-LOC
      real(8), pointer :: Afp(:,:,:)
      real(8), pointer :: Bfp(:,:,:)
      real(8), pointer :: Cfp(:,:,:)
      mwSize, pointer :: dims(:)
      mwSize dimz(3)
      mwSize i
      integer(4) k
      character(len=100) line
!-----
!\
! Check the input
!/
      if( nrhs /= 1 ) then
          call mexErrMsgTxt("Need one 3D double input")
      endif
      if( nlhs > 1 ) then
          call mexErrMsgTxt("Too many outputs")
      endif
!\
! Get a Fortran pointer to the input variable data area.
!/
      Afp => fpGetPr3(prhs(1))
      if( .not.associated(Afp) ) then
          call mexErrMsgTxt("Input is not 3D double")
      endif
!\
! Generate output variable and get pointer to data area
!/
      plhs(1) = mxDuplicateArray(prhs(1))
      Bfp => fpGetPr3(plhs(1))
      if( .not.associated(Bfp) ) then
          call mexErrMsgTxt("Unable to get pointer to output")
      endif
!\
! Perform an array calculation
!/
      k = mexPrintf("Performing OUT = sin(IN) + 7"//achar(10))
      Bfp = sin(Afp) + 7.d0
!\
! Get the dimensions and print them out
!/
      dims => fpGetDimensions(prhs(1))
      k = mexPrintf("Dimensions are:")
      do i=1,size(dims)-1
          write(line,*) dims(i)
          k = mexPrintf(" "//trim(adjustl(line))//" x")
      enddo
      write(line,*) dims(size(dims))
      k = mexPrintf(" "//trim(adjustl(line))//achar(10))
!\
! Allocate memory to make a copy
!/
      dimz = 1
      dimz(1:size(dims)) = dims
      Cfp => fpAllocate(dimz(1),dimz(2),dimz(3))
      if( .not.associated(Cfp) ) then
          call mexErrMsgTxt("Unable to allocate memory for pointer")
      endif
!\
! Make a copy of the input
!/
      Cfp = Afp
!\
! Print out a couple of the rows to demonstrate array slicing
!/
      k = mexPrintf('Print two rows from copy of input'//achar(10))
      write(line,*) 'Copy(1,:,1) = ',Cfp(1,:,1)
      k = mexPrintf(trim(line)//achar(10))
      if( dimz(1) >= 2 .and. dimz(3) >= 2 ) then
          write(line,*) 'Copy(2,:,2) = ',Cfp(2,:,2)
          k = mexPrintf(trim(line)//achar(10))
      endif
!\
! Free the dynamically allocated memory
!/
      call fpDeallocate(Cfp)
!\
! Done
!/
      end subroutine mexFunction