<!DOCTYPE HTML PUBLIC "-//W3C//DTD HTML 4.0 Transitional//EN">
<HTML><HEAD>
<META http-equiv=Content-Type content="text/html; charset=iso-8859-1">
<META content="MSHTML 6.00.2800.1522" name=GENERATOR>
<STYLE></STYLE>
</HEAD>
<BODY bgColor=#ffffff>
<DIV><FONT face=Arial size=2><FONT face="Times New Roman" size=3>Hi
Julien,<BR><BR>Many thanks for your reply. Would you like to have a look at my
requirement?<BR><BR>I'd like to build a spatial relationship between a 3D mesh
(tetrahedron) and<BR>a 3D image. Since the mesh was extracted from CT images
using a software, I<BR>would like to map the mesh to the orignial image data in
my program, so that<BR>I can access all the voxels inside each elements. I don't
know to build a<BR>common physical coordinate system for these two
objects.<BR><BR>The following is my code and running result, the test point [50,
20, 0]<BR>should not in the tetrhedron with<BR>node
no
coordinate<BR>6
65.8246 14.8065
0<BR>1
40.8246 39.8065
25<BR>3
40.8246 14.8065
25<BR>5
65.8246 14.8065 25<BR>however the return bool value of IsInside is
true.<BR><BR>----------------------My<BR>code------------------------------------------------------------------------<BR>---<BR>typedef
itk::DefaultDynamicMeshTraits< float, 3, 3 > MeshTrait;<BR>typedef
itk::Mesh< MeshPixelType, 3, MeshTrait> MeshType;<BR>typedef
MeshType::CellTraits CellTraits;<BR>typedef
itk::CellInterface< float, CellTraits >
CellInterfaceType;<BR>typedef itk::TetrahedronCell< CellInterfaceType
> TetraCellType;<BR>typedef
MeshType::PointType PointType;<BR>typedef
MeshType::CellType CellType;<BR>typedef
CellType::CellAutoPointer
CellAutoPointer;<BR><BR><BR> //access points in cells in a mesh<BR>
typedef MeshType::CellsContainer::Iterator CellIterator;<BR> CellIterator
cellIterator = m_mesh->GetCells()->Begin();<BR> CellIterator end =
m_mesh->GetCells()->End();<BR><BR> // Create a mesh with one
tetrahedtron cell<BR> // in order to extract a bounding box for each
cell<BR> MeshType::Pointer cellMesh = MeshType::New();<BR>
cellMesh->SetCellsAllocationMethod(<BR>
MeshType::CellsAllocatedDynamicallyCellByCell );<BR> CellAutoPointer
oneCell;<BR> oneCell.TakeOwnership( new TetraCellType );<BR>
CellAutoPointer testCell;<BR> testCell.TakeOwnership( new TetraCellType
);<BR><BR> typedef itk::MeshSpatialObject< MeshType >
MeshSpatialObjectType;<BR>
MeshSpatialObjectType::Pointer
cellMeshSpatialObject;<BR><BR> unsigned int numberOfCells =
m_mesh->GetNumberOfCells();<BR><BR> for(unsigned int cellId = 0; cellId
< numberOfCells; cellId++)<BR> {<BR> std::cout << std::endl
<< "Cell ID = " << cellId <<
std::endl;<BR><BR> MeshType::CellType *cellptr1 =
cellIterator.Value();<BR> TetraCellType *cellTetra1 =
dynamic_cast<TetraCellType *>( cellptr1
);<BR> TetraCellType::PointIdIterator pit =
cellTetra1->PointIdsBegin();<BR> TetraCellType::PointIdIterator pitEnd =
cellTetra1->PointIdsEnd();<BR><BR><BR> int id = 0;<BR> while(pit !=
pitEnd)<BR> {<BR> MeshType::PointType pp;<BR>
m_mesh->GetPoint( *pit-1, &pp );<BR> cout << *pit
<< " ";<BR> cout << pp[0] << " " <<
pp[1] << " " << pp[2] << endl;<BR><BR><BR>
//create a mesh with one tetrahedron cell<BR> cellMesh->SetPoint(
*pit-1, pp );<BR> oneCell->SetPointId( id,
*pit);<BR><BR> pit++;<BR>
id++;<BR> }<BR><BR> cellMesh->SetCell( 0, oneCell
);<BR><BR> // print the point info in the mesh with one cell -
tetrahedron<BR> cellMesh->GetCell( 0, testCell
);<BR> TetraCellType::PointIdIterator pit2 =
testCell->PointIdsBegin();<BR> TetraCellType::PointIdIterator pitEnd2 =
testCell->PointIdsEnd();<BR><BR> cout << "cell mesh------"
<< endl;<BR> while(pit2 != pitEnd2)<BR> {<BR>
MeshType::PointType pp1;<BR> cellMesh->GetPoint( *pit2-1,
&pp1 );<BR><BR> cout << *pit2 << "
";<BR> cout << pp1[0] << " " << pp1[1] << "
" << pp1[2] << endl;<BR><BR>
pit2++;<BR> }<BR><BR> // extract a bounding box for the mesh object -
one tetra cell<BR> cellMeshSpatialObject =
MeshSpatialObjectType::New();<BR> cellMeshSpatialObject->SetMesh(
cellMesh );<BR><BR> std::cout << "Cell mesh bounds "
<<<BR> cellMeshSpatialObject->GetBoundingBox()->GetBounds()
<<
std::endl;<BR><BR> //cellMeshSpatialObject->Print(std::cout);<BR><BR><BR> MeshSpatialObjectType::PointType
myPhysicalPoint;<BR> myPhysicalPoint[0] = 50;<BR> myPhysicalPoint[1] =
20;<BR> myPhysicalPoint[2] = 0;<BR><BR> cout << "Is my physical
point inside? : " <<<BR> cellMeshSpatialObject->IsInside(
myPhysicalPoint ) << endl;<BR>
}<BR>----------------------------------------------------------------------------<BR>-----<BR><BR>--------------------------running<BR>result-----------------------------------------<BR>Cell
ID = 0<BR>7 65.8246 39.8065 25<BR>1 40.8246 39.8065 25<BR>6
65.8246 14.8065 0<BR>5 65.8246 14.8065 25<BR>cell mesh------<BR>7
65.8246 39.8065 25<BR>1 40.8246 39.8065 25<BR>6 65.8246 14.8065
0<BR>5 65.8246 14.8065 25<BR>Cell mesh
bounds [40.8246, 65.8246, 14.8065,
39.8065, 0, 25]<BR>Is my physical point [50, 20, 0] inside? : 0<BR><BR>Cell ID =
1<BR>6 65.8246 14.8065 0<BR>3 40.8246 14.8065 25<BR>1 40.8246
39.8065 25<BR>2 40.8246 39.8065 0<BR>cell mesh------<BR>6 65.8246
14.8065 0<BR>3 40.8246 14.8065 25<BR>1 40.8246 39.8065 25<BR>2
40.8246 39.8065 0<BR>Cell mesh bounds
[40.8246, 65.8246, 14.8065, 39.8065, 0, 25]<BR>Is my physical point [50, 20, 0]
inside? : 0<BR><BR>Cell ID = 2<BR>7 65.8246 39.8065 25<BR>6 65.8246
14.8065 0<BR>1 40.8246 39.8065 25<BR>2 40.8246 39.8065 0<BR>cell
mesh------<BR>7 65.8246 39.8065 25<BR>6 65.8246 14.8065 0<BR>1
40.8246 39.8065 25<BR>2 40.8246 39.8065 0<BR>Cell mesh
bounds [40.8246, 65.8246, 14.8065,
39.8065, 0, 25]<BR>Is my physical point [50, 20, 0] inside? : 0<BR><BR>Cell ID =
3<BR>6 65.8246 14.8065 0<BR>2 40.8246 39.8065 0<BR>4 40.8246
14.8065 0<BR>3 40.8246 14.8065 25<BR>cell mesh------<BR>6 65.8246
14.8065 0<BR>2 40.8246 39.8065 0<BR>4 40.8246 14.8065 0<BR>3
40.8246 14.8065 25<BR>Cell mesh bounds
[40.8246, 65.8246, 14.8065, 39.8065, 0, 25]<BR>Is my physical point [50, 20, 0]
inside? : 0<BR><BR>Cell ID = 4<BR>6 65.8246 14.8065 0<BR>1 40.8246
39.8065 25<BR>3 40.8246 14.8065 25<BR>5 65.8246 14.8065 25<BR>cell
mesh------<BR>6 65.8246 14.8065 0<BR>1 40.8246 39.8065 25<BR>3
40.8246 14.8065 25<BR>5 65.8246 14.8065 25<BR>Cell mesh
bounds [40.8246, 65.8246, 14.8065,
39.8065, 0, 25]<BR>Is my physical point [50, 20, 0] inside? : 1<BR><BR>Cell ID =
5<BR>6 65.8246 14.8065 0<BR>2 40.8246 39.8065 0<BR>7 65.8246
39.8065 25<BR>8 65.8246 39.8065 0<BR>cell mesh------<BR>6 65.8246
14.8065 0<BR>2 40.8246 39.8065 0<BR>7 65.8246 39.8065 25<BR>8
65.8246 39.8065 0<BR>Cell mesh bounds
[40.8246, 65.8246, 14.8065, 39.8065, 0, 25]<BR>Is my physical point [50, 20, 0]
inside? :
0<BR>----------------------------------------------------------------------------<BR>---------------------<BR><BR>Thanks
a lot.<BR><BR>Kind regards,<BR><BR>Yan<BR><BR><BR>----- Original Message -----
<BR>From: "Julien Jomier" <</FONT><A href="mailto:jjomier@cs.unc.edu"><FONT
face="Times New Roman" size=3>jjomier@cs.unc.edu</FONT></A><FONT
face="Times New Roman" size=3>><BR>To: <</FONT><A
href="mailto:y.yang@anglia.ac.uk"><FONT face="Times New Roman"
size=3>y.yang@anglia.ac.uk</FONT></A><FONT face="Times New Roman"
size=3>><BR>Sent: Friday, December 09, 2005 6:29 PM<BR>Subject: Re:
IsInside() method in itk::MeshSpatialObject<BR><BR><BR>> Hi
Yan,<BR>><BR>> Can you be more specific? and send me a sample code that
shows the<BR>problem?<BR>> I cannot reproduce the problem you are
describing.<BR>><BR>> Let me know,<BR>><BR>> Julien<BR>><BR>>
Yan Yang wrote:<BR>> > Hi Julien,<BR>> ><BR>> > I am using
ITK2.4.1. I am trying to query if a point is inside a<BR>> > tetrahedron.
I created a Mesh Spatial Object with only one tetrahedron<BR>cell,<BR>> >
and used IsInside() method to query if a point within the
cell's<BR>bounding<BR>> > box is inside the tetrahedron. But I am afraid
there might be something<BR>> > wrong with the return boolean value. I
would be grateful if anyone could<BR>> > tell me if the IsInside()
function can check a point inside the mesh.<BR>> ><BR>> > Kind
regards,<BR>> ><BR>> > Yan<BR>> ><BR>> ><BR>>
></FONT><BR></FONT></DIV></BODY></HTML>