"-Ob","-O[utput]b[lack]","Output: black","Optional","1","1","OutBlack","outFile::mrcImage","NULL"
"-Owe","-O[utput]w[hite]e[dge]","Output: whiteEdge","Optional","1","1","OutWhiteEdge","outFile::mrcImage","NULL"
"-Obe","-O[utput]b[lack]e[dge]","Output: blackEdge","Optional","1","1","OutBlackEdge","outFile::mrcImage","NULL"
+"-Odir","-O[utput]dir[ection]","Output: basename.[0|1|2] vector for any","Optional","1","1","Direction","String","direction"
"-lt","-l[ine]t[hickness]","[pixel]","Optional","2","1","LineThicknessX","Real","1","2","LineThicknessY","Real","1"
"-lm","-l[ine]m[ode]","LineMode","Optional","1","1","LineMode","Integer","0"
"-ls","-l[ine]s[hape]","LineShape","Optional","1","1","LineShape","Integer","0"
+"-deltaPhi","-deltaPhi","deltaPhi","Optional","1","1","deltaPhi","Real","30"
+"-deltaTheta","-deltaTheta","deltaTheta","Optional","1","1","deltaTheta","Real","30"
"-c","-c[onfig]","ConfigurationFile","Optional","1","1","configFile","inFile","NULL"
"-m","-m[ode]","Mode","Optional","1","1","mode","Integer","0"
mrcImage black;
mrcImage whiteEdge;
mrcImage blackEdge;
+ //
+ mrcImage direction[3];
+
// Control: Input
mrcImageParaTypeRealCoord Length;
mrcImageParaTypeRealCoord LineThickness;
lmrcImageMeanFreePathCalcLineMode LineMode;
lmrcImageMeanFreePathCalcLineShape LineShape;
+ mrcImageParaTypeReal deltaPhi;
+ mrcImageParaTypeReal deltaTheta;
// TemporalVariables
mrcImageParaTypeRealCoord min;
MeanFreePathCalcModeAlongZaxis=0,
MeanFreePathCalcModeAlongXaxis=1,
MeanFreePathCalcModeAlongYaxis=2,
+ MeanFreePathCalcModeAnyDirection=3,
MeanFreePathCalcModeCalcMode=0x0f,
MeanFreePathCalcModeRealLength=0x10
} lmrcImageMeanFreePathCalcMode;
extern void lmrcImageMeanFreePathCalcPrint(FILE* fpt, lmrcImageMeanFreePathCalcInfo* linfo, lmrcImageMeanFreePathCalcMode mode);
extern void lmrcImageMeanFreePathCalc(mrcImage* in, lmrcImageMeanFreePathCalcInfo* linfo, lmrcImageMeanFreePathCalcMode mode);
extern void lmrcImageMeanFreePathCalcAlongZaxis(mrcImage* in, lmrcImageMeanFreePathCalcInfo* linfo, lmrcImageMeanFreePathCalcMode mode);
+extern void lmrcImageMeanFreePathCalcAnyDirection(mrcImage* in, lmrcImageMeanFreePathCalcInfo* linfo, lmrcImageMeanFreePathCalcMode mode);
extern double lmrcImageMeanFreePathCalcBWEvaluation(mrcImage* in, lmrcImageMeanFreePathCalcInfo* linfo, lmrcImageMeanFreePathCalcMode mode);
+extern void lmrcImageMeanFreePathCalcMeanCalc(mrcImage* in, lmrcImageMeanFreePathCalcInfo* linfo, lmrcImageMeanFreePathCalcMode mode);
extern void lmrcImageMeanFreePathCalcModePrint(FILE* fpt);
mrcImageMeanFreePathCalcInfo info;
lmrcImageMeanFreePathCalcInfo linfo;
mrcImage in;
+ int i;
init0(&info);
argCheck(&info, argc, argv);
linfo.LineMode = info.LineMode;
linfo.LineShape = info.LineShape;
+ linfo.deltaPhi = info.deltaPhi*RADIAN;
+ linfo.deltaTheta = info.deltaTheta*RADIAN;
+
DEBUGPRINT1("mode: %ld\n", info.mode);
lmrcImageMeanFreePathCalc(&in, &linfo, info.mode);
if(info.flagOutBlackEdge) {
mrcFileWrite(&linfo.blackEdge, info.OutBlackEdge, "in main", 0);
}
+
+ if(info.flagDirection && info.mode == 3) {
+ char s[1024];
+ for(i=0; i<3; i++) {
+ sprintf(s, "%s.%d", info.Direction, i);
+ mrcFileWrite(&linfo.direction[i], s, "in main", 0);
+ }
+ }
return EXIT_SUCCESS;
}
lmrcImageMeanFreePathCalcModePrint(FILE* fpt)
{
fprintf(fpt, ">>>> mode \n");
- fprintf(fpt, "%d: AlongZ-axis\n", MeanFreePathCalcModeAlongZaxis);
- fprintf(fpt, "%d: AlongX-axis\n", MeanFreePathCalcModeAlongXaxis);
- fprintf(fpt, "%d: AlongY-axis\n", MeanFreePathCalcModeAlongYaxis);
+ fprintf(fpt, "%d: Along Z-axis\n", MeanFreePathCalcModeAlongZaxis);
+ fprintf(fpt, "%d: Along X-axis\n", MeanFreePathCalcModeAlongXaxis);
+ fprintf(fpt, "%d: Along Y-axis\n", MeanFreePathCalcModeAlongYaxis);
+ fprintf(fpt, "%d: For any direction\n", MeanFreePathCalcModeAnyDirection);
fprintf(fpt, "+%d: Length is shown in real using HeaderLength \n", MeanFreePathCalcModeRealLength);
fprintf(fpt, ">>>> line shape\n");
lmrcImageMeanFreePathCalcAlongZaxis(&tmp, linfo, mode);
break;
}
+ case MeanFreePathCalcModeAnyDirection: {
+ DEBUGPRINT1("mode: %d in Any Direction\n", mode);
+ lmrcImageMeanFreePathCalcAnyDirection(in, linfo, mode);
+ break;
+ }
default: {
DEBUGPRINT1("mode: %d in default\n", mode);
fprintf(stderr, "Not supported mode: %d\n", mode);
break;
}
}
-
}
void
}
}
}
+ lmrcImageMeanFreePathCalcMeanCalc(in, linfo, mode);
+}
+
+void
+lmrcImageMeanFreePathCalcMeanCalc(mrcImage* in, lmrcImageMeanFreePathCalcInfo* linfo, lmrcImageMeanFreePathCalcMode mode)
+{
+ int i;
linfo->averageWhite = linfo->averageBlack = 0;
linfo->averageWhiteEdge = linfo->averageBlackEdge = 0;
return data;
}
+
+void
+lmrcImageMeanFreePathCalcAnyDirection(mrcImage* in, lmrcImageMeanFreePathCalcInfo* linfo, lmrcImageMeanFreePathCalcMode mode)
+{
+ /*
+ mrcImageParaTypeRealCoord Length;
+ mrcImageParaTypeRealCoord LineThickness;
+ lmrcImageMeanFreePathCalcLineMode LineMode;
+ lmrcImageMeanFreePathCalcLineShape LineShape;
+ */
+
+ mrcImageParaTypeReal x, y, z;
+ mrcImageParaTypeReal dstx, dsty, dstz;
+ mrcImageParaTypeReal phi, theta, dphi, dtheta, dr, r, rmax, rmaxmax, thetamax, phimax;
+ int ntheta, nphi, itheta, iphi, ir, nr;
+ int flagEdge, flagEdgeMax;
+ double data, dstdata;
+ int i, irmaxmax;
+ int* countWhite;
+ int* countBlack;
+ int* countWhiteEdge;
+ int* countBlackEdge;
+
+ linfo->white.Header = in->Header;
+ linfo->black.Header = in->Header;
+ linfo->whiteEdge.Header = in->Header;
+ linfo->blackEdge.Header = in->Header;
+ mrcInit(&linfo->white, NULL);
+ mrcInit(&linfo->black, NULL);
+ mrcInit(&linfo->whiteEdge, NULL);
+ mrcInit(&linfo->blackEdge, NULL);
+
+ for(i=0; i<3; i++) {
+ linfo->direction[i].Header = in->Header;
+ mrcInit(&linfo->direction[i], NULL);
+ }
+
+ nphi = floor(2*M_PI/linfo->deltaPhi + 1);
+ if(1<in->HeaderN.z) {
+ ntheta = floor(M_PI/linfo->deltaTheta + 1);
+ } else {
+ ntheta = 1;
+ theta = 0;
+ }
+ nr = sqrt(SQR(in->HeaderN.x)+SQR(in->HeaderN.y)+SQR(in->HeaderN.z));
+
+ DEBUGPRINT3("(nr, nphi, ntheta) = (%d %d %d)\n", nr, nphi, ntheta);
+ countBlack = memoryAllocate(sizeof(int)*(nr+1), "in lmrcImageMeanFreePathCalc");
+ countWhite = memoryAllocate(sizeof(int)*(nr+1), "in lmrcImageMeanFreePathCalc");
+ countBlackEdge = memoryAllocate(sizeof(int)*(nr+1), "in lmrcImageMeanFreePathCalc");
+ countWhiteEdge = memoryAllocate(sizeof(int)*(nr+1), "in lmrcImageMeanFreePathCalc");
+
+ linfo->countBlack = countBlack;
+ linfo->countWhite = countWhite;
+ linfo->countBlackEdge = countBlackEdge;
+ linfo->countWhiteEdge = countWhiteEdge;
+ linfo->N = nr;
+
+ for(z=0; z<in->HeaderN.z; z++) {
+ for(y=0; y<in->HeaderN.y; y++) {
+ for(x=0; x<in->HeaderN.x; x++) {
+
+ mrcPixelDataGet(in, x, y, z, &data, mrcPixelRePart, mrcPixelHowNearest);
+
+ rmaxmax = 0;
+ for(itheta=0; itheta<ntheta; itheta++) {
+ if(1<in->HeaderN.z) {
+ theta = itheta*linfo->deltaTheta - M_PI/2.0;
+ } else {
+ theta = 0;
+ }
+ for(iphi=0; iphi<nphi; iphi++) {
+ phi = iphi*linfo->deltaPhi;
+
+ flagEdge = 0;
+ rmax = -1;
+ for(ir=1; ir<nr; ir++) {
+ r=ir;
+
+ dstx = r*cos(phi)*cos(theta)+x;
+ dsty = r*sin(phi)*cos(theta)+y;
+ dstz = r*sin(theta) +z;
+
+ if(-0.5<=dstx && dstx < in->HeaderN.x -0.5
+ &&-0.5<=dsty && dsty < in->HeaderN.y -0.5
+ &&-0.5<=dstz && dstz < in->HeaderN.z -0.5) {
+ mrcPixelDataGet(in, dstx, dsty, dstz, &dstdata, mrcPixelRePart, mrcPixelHowNearest);
+ } else {
+ rmax = r - 1;
+ flagEdge = 1;
+ break;
+ }
+ if(dstdata != data) {
+ rmax = r - 1;
+ break;
+ }
+ }
+ if(rmaxmax<rmax) {
+ rmaxmax = rmax;
+ phimax = phi;
+ thetamax = theta;
+ flagEdgeMax = flagEdge;
+ }
+ }
+ }
+
+ irmaxmax = (int)rmaxmax;
+ if(0<data) { // white
+ if(0<flagEdgeMax) { // whiteEdge
+ mrcPixelDataSet(&linfo->whiteEdge, x, y, z, rmaxmax, mrcPixelRePart);
+ countWhiteEdge[irmaxmax]++;
+ } else { // white
+ mrcPixelDataSet(&linfo->white, x, y, z, rmaxmax, mrcPixelRePart);
+ countWhite[irmaxmax]++;
+ }
+ } else { // black
+ if(0<flagEdgeMax) { // blackEdge
+ mrcPixelDataSet(&linfo->blackEdge, x, y, z, rmaxmax, mrcPixelRePart);
+ countBlackEdge[irmaxmax]++;
+ } else { // black
+ mrcPixelDataSet(&linfo->black, x, y, z, rmaxmax, mrcPixelRePart);
+ countBlack[irmaxmax]++;
+ }
+ }
+ mrcPixelDataSet(&linfo->direction[0], x, y, z, cos(phimax)*cos(thetamax), mrcPixelRePart);
+ mrcPixelDataSet(&linfo->direction[1], x, y, z, sin(phimax)*cos(thetamax), mrcPixelRePart);
+ mrcPixelDataSet(&linfo->direction[2], x, y, z, sin(thetamax), mrcPixelRePart);
+ }
+ }
+ }
+
+ lmrcImageMeanFreePathCalcMeanCalc(in, linfo, mode);
+}