New Fortran STL Binary File Reader

The old one still wins


Well Dan won’t be happy with it. I am not very happy with it. But While thinking about optimising the old fortran binary stl reader for speed, I knew the speed can be reduced significantly if I can reduce the number of I/O calls in the routine. For that need a mechanism to get similar functionality like seek and tell in fortran.

After researching it, found this good notes on the new stream I/O feature of modern fortran. And using Dr. Clive Page’s some notes on Stream I/O., built this new version of binary stl file reader. The present code is compact, straight-forward and easy to understand but isn’t speedier.

so the Quest continues.


       program stlread
      
       character header*80,filname*80
       real*4 n(3),x1(3),x2(3),x3(3)
       integer*2 padding
       integer*4 ntri,iunit
       real*4, allocatable :: normals(:,:),triangles(:,:)
       integer iargc

       if (iargc() .ge. 1) then
           call getarg(1,filname)
       else
           filname="/Users/Sukhbinder/Desktop/projects/projects/work/stlreader/zDino/Dino.stl"
       end if

       iunit=13
       print *, filname
       open(unit=iunit,file=filname,status='old',access='stream', &
     & form='unformatted')

       read(iunit)header
       read(iunit) ntri
       
       allocate(normals(3,ntri))
       allocate(triangles(3,ntri*3))

       k=1
       do i = 1,ntri
        read(iunit) normals(1,i),normals(2,i),normals(3,i)
        read(iunit) triangles(1,k),triangles(2,k),triangles(3,k)
        read(iunit) triangles(1,k+1),triangles(2,k+1),triangles(3,k+1)
        read(iunit) triangles(1,k+2),triangles(2,k+2),triangles(3,k+2)
        read(iunit) padding
        k=k+3
       end do
       close(iunit)
       print *, header
       print *, ntri
       print *, normals(:,1)
       print *, triangles(:,1)
       print *, triangles(:,2)
       print *, triangles(:,3)

       end 

And here’s the old one.


    program stlread

! Program to Read Binary STL file.
! programmed by Sukhbinder Singh
! stores normals and the triangles in the normals and triangles arrays.
!
! 23rd May 2011

      real*4 n(3),x1(3),x2(3),x3(3)
      integer*2 padding
      integer*4 ntri,iunit,irc
      character(len=80) title, filename
      character(len=52) sentry
      double precision, allocatable :: normals(:,:),triangles(:,:)
      integer iargc

      if (iargc() .ge. 1) then
        call getarg(1,filename)
      else
        filename="/Users/Sukhbinder/Desktop/projects/projects/work/stlreader/zDino/Dino.stl"
      end if

      iunit=13
     
      open(unit=iunit,file=filename,status='old',form='unformatted', &
    & access='direct',recl=80)

      read(iunit,rec=1)title
      close(iunit,status='keep')

      open(unit=iunit,file=filename,status='old',form='unformatted', &
    & access='direct',recl=4)
     
      read(iunit,rec=21)ntri
      close(iunit,status='keep')
           
      open(unit=iunit,file=filename,status='old',form='unformatted', &
    & access='direct',recl=2)
      irc=(80+4)/2+1
     
      !call fseek(iunit,84,0,ierr)
     
      allocate(normals(3,ntri))
      allocate(triangles(3,ntri*3))
      
      !read(iunit,rec=1)sentry
      !print *, sentry
      k=1
      do i=1,2 !ntri
       call readbin(iunit,n(1),irc) 
       call readbin(iunit,n(2),irc)
       call readbin(iunit,n(3),irc)

       call readbin(iunit,x1(1),irc)
       call readbin(iunit,x1(2),irc)
       call readbin(iunit,x1(2),irc)

       call readbin(iunit,x2(1),irc)
       call readbin(iunit,x2(2),irc)
       call readbin(iunit,x2(3),irc)

       call readbin(iunit,x3(1),irc)
       call readbin(iunit,x3(2),irc)
       call readbin(iunit,x3(3),irc)

       normals(1,i)=n(1)
       normals(2,i)=n(2)
       normals(3,i)=n(3)
     

       triangles(1,k)=x1(1)
       triangles(2,k)=x1(2)
       triangles(3,k)=x1(3)

       triangles(1,k+1)=x2(1)
       triangles(2,k+1)=x2(2)
       triangles(3,k+1)=x2(3)


       triangles(1,k+2)=x3(1)
       triangles(2,k+2)=x3(2)
       triangles(3,k+2)=x3(3)

       if (i .eq. 1) then
         write(*,*) "normals: ",n
         write(*,*) "v1 : ",x1
         write(*,*) "v2 : ",x2
         write(*,*) "v3 : ",x3

       end if
       k=k+3
       read(iunit,rec=irc)padding
       irc=irc+1
      end do   

      write(*,*) trim(filename),' has this title ',trim(title),' and has ',ntri, ' triangles' 
      end program

      subroutine readbin(iunit,a,irc)
      real*4 a,b
      integer*2 ib(2)
      equivalence(b,ib)

      read(iunit,rec=irc)ib(1)
      irc=irc+1
      read(iunit,rec=irc)ib(2)
      irc=irc+1
      
      a=b

      end

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out / Change )

Twitter picture

You are commenting using your Twitter account. Log Out / Change )

Facebook photo

You are commenting using your Facebook account. Log Out / Change )

Google+ photo

You are commenting using your Google+ account. Log Out / Change )

Connecting to %s