Changes between Version 4 and Version 5 of CameraObjects

Timestamp:

04/26/12 10:42:12 (13 years ago)

Author:

Jonathan

Comment:

—

CameraObjects

@@ -79 +79 @@
    END SUBROUTINE BinCell
 }}}
- We could in principal continue to raise the sample_res but it gets computationally expensive... A better approach would be to sample along each ray at some interval...  Then the sampling is pixel based instead of volume based - and should reduce the aliasing...
+ We could in principal continue to raise the sample_res but it gets computationally expensive... A better approach would be to sample along each ray at some interval...  Then the sampling is pixel based instead of volume based - and should reduce the aliasing...  Or we could calculate the optical depth across each cell for each ray that intersects the cell.  Given that our mesh is cartesian this is not as hard as it might sound.
+
+ * Find what rays could intersect the cell
+ * For each ray find what which two faces the ray intersects and the corresponding points of entry and exit
+ * The contribution to the pixel will be the the value of the cell times the distance between those points.
+
+{{{
+
+   SUBROUTINE BinCell(Camera, pos, dx, data, rho)
+      TYPE(CameraDef), POINTER :: Camera
+      REAL(KIND=qPREC) :: pos(3),xpos(3)
+      REAL(KIND=qPREC) :: dx, ddx, xbounds(3,2), my_pixel(2)
+      REAL(KIND=qPREC), DIMENSION(:,:) :: data
+      REAL(KIND=qPREC) :: rho,a,intersection(6,3), ray(3), max_distance
+      INTEGER :: ipos(2), sample_res, i, j, k, npoints, min_pixels(2), max_pixels(2), pixel(2), dim, odim(2), edge
+      xbounds(:,1)=(pos-half*dx)
+      xbounds(:,2)=(pos+half*dx)
+      min_pixels=huge(min_pixels(1))
+      max_pixels=0
+      DO i=1,2
+         DO j=1,2
+            DO k=1,2
+               my_pixel=GetPixel(Camera, (/xbounds(1,i), xbounds(2,j), xbounds(3,k)/))*shape(data)+half
+               min_pixels=max(1,min(min_pixels, floor(my_pixel)))
+               max_pixels=min(shape(data),max(max_pixels, ceiling(my_pixel)))
+            END DO
+         END DO
+      END DO
+      DO i=min_pixels(1), max_pixels(1)
+         DO j=min_pixels(2), max_pixels(2)
+            pixel=(/i,j/)-half
+            Ray=GetRay(Camera, REAL(pixel,KIND=qPREC)/REAL(shape(data), KIND=qPREC))
+            npoints=0
+            DO dim=1,3
+               DO edge=1,2
+                  ! Camera%pos(dim)+a*ray(dim)=xbounds(dim,edge)
+                  a=(xbounds(dim,edge)-Camera%pos(dim))/ray(dim)
+                  xpos=Camera%pos+a*ray
+                  odim=modulo((/dim,dim+1/),3)+1
+                  IF (ALL(xpos(odim) >= xbounds(odim,1) .AND. xpos(odim) <= xbounds(odim,2))) THEN
+                     npoints=npoints+1
+                     intersection(npoints,:)=xpos
+                  END IF
+               END DO
+            END DO
+            IF (npoints == 0) CYCLE
+            max_distance=0d0
+            DO k=1,npoints
+               max_distance=max(max_distance, sqrt(sum((intersection(k,:)-intersection(1,:))**2)))
+            END DO
+            data(i,j)=data(i,j)+rho*max_distance
+         END DO
+      END DO
+   END SUBROUTINE BinCell
+}}}