// -*- FreeFem++ -*- // Time-stamp: "2020-08-04 17:56:11 fujiwara" // // (Original) // Larry Shepp, Benjamin F. Logan, // The Fourier reconstruction of a head section, // IEEE Trans. Nucl. Sci. Vol. NS-21(3) pp.21--43 (June 1974). // // Modification (in contrast) // Toft, Peter Aundal (Toft, P.A.), // The Radon Transform - Theroy and Implementation, // Ph.D. Thesis, Department of Mathematical Modelling, // Technical University of Denmark (1996) // (Table B.1 p.201) // //---------------------------------------------------------------------- // Domain //---------------------------------------------------------------------- //border RB(t=-1,1) { x=t; y=-1; label=1; }; //border RR(t=-1,1) { x=1; y=t; label=2; }; //border RT(t=1,-1) { x=t; y=1; label=3; }; //border RL(t=1,-1) { x=-1; y=t; label=4; }; border a(t=0,2*pi) { x=0.69*cos(t); y=0.92*sin(t); label=1; }; border b(t=0,2*pi) { x=0.6624*cos(t); y=0.874*sin(t)-0.0184; label=2; }; real tt=18.0/360.0*2*pi; real c1=cos(-tt); real s1=sin(-tt); border c(t=0,2*pi) { x= c1*(0.11*cos(t)) - s1*(0.31*sin(t))+0.22; y= s1*(0.11*cos(t)) + c1*(0.31*sin(t)); label=3; }; real c2=cos(tt); real s2=sin(tt); border d1(t=5.973076536681481-2*pi,0.206321103590681) { x = c2*(0.16*cos(t)) - s2*(0.41*sin(t)) - 0.22; y = s2*(0.16*cos(t)) + c2*(0.41*sin(t)); label=4; }; border d2(t=0.206321103590681,0.8173443749584077) { x = c2*(0.16*cos(t)) - s2*(0.41*sin(t)) - 0.22; y = s2*(0.16*cos(t)) + c2*(0.41*sin(t)); label=4; }; border d3(t=0.8173443749584077,5.799123061766094) { x = c2*(0.16*cos(t)) - s2*(0.41*sin(t)) - 0.22; y = s2*(0.16*cos(t)) + c2*(0.41*sin(t)); label=4; }; border d4(t=5.799123061766094,5.973076536681481) { x = c2*(0.16*cos(t)) - s2*(0.41*sin(t)) - 0.22; y = s2*(0.16*cos(t)) + c2*(0.41*sin(t)); label=4; }; border e1(t=4.931330217751022-2*pi,3.269068382966671) { x=0.21*cos(t); y=0.25*sin(t)+0.35; label=5; }; border e2(t=3.269068382966671, 4.232175040320118) { x=0.21*cos(t); y=0.25*sin(t)+0.35; label=5; }; border e3(t=4.232175040320118, 4.493447743018358) { x=0.21*cos(t); y=0.25*sin(t)+0.35; label=5; }; border e4(t=4.493447743018358, 4.931330217751022) { x=0.21*cos(t); y=0.25*sin(t)+0.35; label=5; }; border f1(t=3.011486528670876-2*pi,0.130106124918917) { x=0.046*cos(t); y=0.046*sin(t)+0.1; label=6; }; border f2(t=0.130106124918917, 3.011486528670876) { x=0.046*cos(t); y=0.046*sin(t)+0.1; label=6; }; border g1(t=4.102403867569372-2*pi,2.484669020407591) { x=0.046*cos(t); y=0.046*sin(t)-0.1; label=7; }; border g2(t=2.484669020407591, 4.102403867569372) { x=0.046*cos(t); y=0.046*sin(t)-0.1; label=7; }; border h(t=0,2*pi) { x=0.046*cos(t)-0.08; y=0.023*sin(t)-0.605; label=8; }; border i(t=0,2*pi) { x=0.023*cos(t); y=0.023*sin(t)-0.605; label=9; }; border j(t=0,2*pi) { x=0.023*cos(t)+0.06; y=0.046*sin(t)-0.605; label=10; }; //---------------------------------------- // length and ratio of each segment //---------------------------------------- // RB,RR,RT,RL : 2.0000 35.5619 // a : 5.0838 90.3950 // b : 4.8497 86.2316 // c : 1.3954 24.8115 // d1 : 0.2095 3.7254 // d2 : 0.2210 3.9295 // d3 : 1.3808 24.5522 // d4 : 0.0666 1.1838 // e1 : 1.0716 19.0542 // e2 : 0.2278 4.0497 // e3 : 0.0562 1.0000 <-- minimum, ratio=1 // e4 : 0.0923 1.6404 // f1 : 0.1565 2.7824 // f2 : 0.1325 2.3568 // g1 : 0.2146 3.8160 // g2 : 0.0744 1.3232 // h : 0.2228 3.9622 // i : 0.1445 2.5696 // j : 0.2228 3.9622 //---------------------------------------- int m = 2; //mesh Th = buildmesh( RB(m)+RR(m)+RT(m)+RL(m)+a(5*m)+b(5*m)+c(5*m)+d1(2*m)+d2(m)+d3(4*m)+d4(m)+e1(3*m)+e2(m)+e3(m)+e4(m)+f1(m)+f2(m)+g1(m)+g2(m)+h(m)+i(m)+j(m) ); mesh Th = buildmesh( a(90*m)+b(86*m) +c(25*m) + d1(4*m)+d2(4*m)+d3(25*m)+d4(m) + e1(19*m)+e2(4*m)+e3(m)+e4(2*m) +f1(3*m)+f2(2*m) + g1(4*m)+g2(m) +h(4*m)+i(3*m)+j(4*m) ); savemesh( Th, "SheppLogan" + m + ".msh" ); load "plotPDF" // without semicolon plotPDF( Th, "SheppLogan"+m, belabel=true, index=false, withmesh=0.1 ); //---------------------------------------------------------------------- // End of file //----------------------------------------------------------------------