diff --git a/Makefile b/Makefile index 8cd90ede..4f264a7b 100644 --- a/Makefile +++ b/Makefile @@ -48,7 +48,7 @@ SOBJ = network.o ifile.o conman.o SCinter.o splitter.o passwd.o \ circular.o el755driv.o maximize.o sicscron.o tecsdriv.o sanscook.o \ tasinit.o tasutil.o t_rlp.o t_conv.o d_sign.o d_mod.o \ tasdrive.o tasscan.o synchronize.o definealias.o swmotor.o t_update.o \ - hmcontrol.o userscan.o slsmagnet.o rs232controller.o + hmcontrol.o userscan.o slsmagnet.o rs232controller.o lomax.o MOTOROBJ = motor.o el734driv.o simdriv.o el734dc.o pipiezo.o pimotor.o COUNTEROBJ = countdriv.o simcter.o counter.o diff --git a/danu.dat b/danu.dat index 265dcb48..12e2756e 100644 --- a/danu.dat +++ b/danu.dat @@ -1,3 +1,3 @@ - 103 + 123 NEVER, EVER modify or delete this file You'll risk eternal damnation and a reincarnation as a cockroach!|n \ No newline at end of file diff --git a/lomax.c b/lomax.c new file mode 100644 index 00000000..99612754 --- /dev/null +++ b/lomax.c @@ -0,0 +1,476 @@ +/*------------------------------------------------------------------------- + L o c a l M a x i m u m + + This is a module for searching a local maximum in a 2D histogram as + collected at TRICS. + + copyright: see copyright.h + + Mark Koennecke, November 2001 +-------------------------------------------------------------------------*/ +#include +#include +#include +#include +#include "fortify.h" +#include "sics.h" +#include "lomax.h" +#include "lomax.i" +#include "HistMem.h" + +#define WINDOW 0 +#define THRESHOLD 1 +#define STEEPNESS 2 +#define COGWINDOW 3 +#define COGCONTOUR 4 + +/*-------------------------------------------------------------------*/ +static int testBoundaries(int xsize, int ysize, int window, + int i, int j) +{ + int half; + + /* + if the window touches the data boundary the result is probably not + very useful. Discarding these cases early on also protects us + against data array overrun problems. + */ + + half = window/2; + + if( i - half < 0 || i + half > xsize) + return 0; + + if( j - half < 0 || j + half > ysize) + return 0; + + return 1; +} +/*-------------------------------------------------------------------*/ +static int testSteepness(int *iData, int xsize, int ysize, + int i, int j, int window, int steepness) +{ + int testValue, x, y, half; + int *iPtr; + + + testValue = iData[j * xsize + i] - steepness; + half = window/2; + + /* + test upper row + */ + iPtr = iData + (j - half) * xsize + i - half; + for(x = 0; x < window; x++) + { + if(iPtr[x] > testValue) + return 0; + } + + /* + test lower row + */ + iPtr = iData + (j + half) * xsize + i - half; + for(x = 0; x < window; x++) + { + if(iPtr[x] > testValue) + return 0; + } + + /* + test columns + */ + for(y = j - half; y < j + half; y++) + { + /* + left + */ + if(iData[y*xsize + i - half] > testValue) + return 0; + /* + right + */ + if(iData[y*xsize + i + half] > testValue) + return 0; + } + + return 1; +} + +/*--------------------------------------------------------------------*/ +static int testMaximum(int *iData, int xsize, int ysize, + int i, int j, int window) +{ + int testValue, x, y, half; + int *iPtr; + + testValue = iData[j * xsize + i]; + half = window/2; + + for(y = j - half; y < j + half; y++) + { + iPtr = iData + y * xsize + i - half; + for(x = 0; x < window; x++) + { + if(iPtr[x] > testValue) + return 0; + } + } + return 1; +} +/*--------------------------------------------------------------------*/ +int testLocalMaximum(int *iData, int xsize, int ysize, + int i, int j, + int window, int steepness, int threshold, + int *intensity) +{ + + if(!testBoundaries(xsize,ysize,window,i,j)) + return 0; + + + if(!testMaximum(iData,xsize,ysize,i,j,window)) + return 0; + + if(iData[j * xsize + i] < threshold) + return 0; + + if(!testSteepness(iData,xsize,ysize,i,j,window, steepness)) + return 0; + + *intensity = iData[j * xsize + i]; + + return 1; +} +/*-------------------------------------------------------------------*/ +int calculateCOG(int *iData, int xsize, int ysize, + int *i, int *j, int *intensity, + int cogWindow, + float contour) +{ + int x, xLow, xMax, y, yLow, yMax; + int threshold; + float cogTotal, cogX, cogY; + int *iPtr; + + /* + preparations + */ + xLow = *i - cogWindow/2; + if(xLow < 0) + xLow = 0; + xMax = *i + cogWindow/2; + if(xLow >= xsize) + xMax = xsize-1; + + yLow = *j - cogWindow/2; + if(yLow < 0) + yLow = 0; + yMax = *j + cogWindow/2; + if(yLow >= ysize) + yMax = ysize-1; + + threshold = (int)(float)iData[*j * xsize + *i] * contour; + + /* + build the sums + */ + cogTotal = cogY = cogX = .0; + for(y = yLow; y < yMax; y++) + { + iPtr = iData + y *xsize; + for(x = xLow; x < xMax; x++) + { + if(iPtr[x] > threshold) + { + cogTotal += iPtr[x]; + cogY += y * iPtr[x]; + cogX += x * iPtr[x]; + } + } + } + if(cogTotal < .0) + { + return 0; + } + + *i = (int)nintf(cogX/cogTotal); + *j = (int)nintf(cogY/cogTotal); + *intensity = (int)cogTotal; + + return 1; +} +/*-------------------------------------------------------------------*/ +static int checkHM(pHistMem *pHM, SicsInterp *pSics, SConnection *pCon, + char *name, int *iDim) +{ + CommandList *pCom = NULL; + char pBueffel[256]; + int nDim; + + pCom = FindCommand(pSics, name); + if(!pCom) + { + sprintf(pBueffel,"ERROR: histogram memory %s not found", name); + SCWrite(pCon,pBueffel,eError); + return 0; + } + if(!pCom->pData) + { + sprintf(pBueffel,"ERROR: histogram memory %s not found", name); + SCWrite(pCon,pBueffel,eError); + return 0; + } + *pHM = (pHistMem)pCom->pData; + if(!iHasType(*pHM, "HistMem")) + { + sprintf(pBueffel,"ERROR: %s is no histogram memory!", name); + SCWrite(pCon,pBueffel,eError); + return 0; + } + /* + we now know that we have a histogram memory, now check Sizes + */ + GetHistDim(*pHM,iDim,&nDim); + if(nDim < 2) + { + sprintf(pBueffel,"ERROR: %s is not 2 dimensional!", name); + SCWrite(pCon,pBueffel,eError); + return 0; + } + + return 1; +} +/*--------------------------------------------------------------------*/ +int LoMaxAction(SConnection *pCon, SicsInterp *pSics, + void *pData, int argc, char *argv[]) +{ + pLoMax self = NULL; + char pBueffel[256], pNum[20]; + int iRet, i, j, intensity; + int iDim[10], nDim; + pHistMem pHM = NULL; + CommandList *pCom = NULL; + Tcl_DString result; + int window; + int *iData; + double dVal; + ObPar *ob = NULL; + + self = (pLoMax)pData; + assert(pCon); + assert(pSics); + assert(self); + + /* + we need arguments + */ + if(argc < 2) + { + sprintf(pBueffel,"ERROR: insufficient number of arguments to %s", + argv[0]); + SCWrite(pCon,pBueffel,eError); + return 0; + } + + /* + interpret arguments + */ + strtolower(argv[1]); + if(strcmp(argv[1],"search") == 0) + { + if(argc < 3) + { + sprintf(pBueffel,"ERROR: insufficient number of arguments to %s.search", + argv[0]); + SCWrite(pCon,pBueffel,eError); + return 0; + } + if(!checkHM(&pHM, pSics,pCon, argv[2],iDim)) + { + return 0; + } + Tcl_DStringInit(&result); + Tcl_DStringStartSublist(&result); + iData = GetHistogramPointer(pHM,pCon); + window = (int)ObVal(self->pParam,WINDOW); + for(i = 0 + window/2; i < iDim[0] - window/2; i++) + { + for(j = 0 + window/2; j < iDim[1] - window/2; j++) + { + if(testLocalMaximum(iData,iDim[0], iDim[1], + i,j, + (int)ObVal(self->pParam,WINDOW), + (int)ObVal(self->pParam,STEEPNESS), + (int)ObVal(self->pParam,THRESHOLD), + &intensity)) + { + Tcl_DStringStartSublist(&result); + sprintf(pNum,"%d", i); + Tcl_DStringAppendElement(&result,pNum); + sprintf(pNum,"%d", j); + Tcl_DStringAppendElement(&result,pNum); + sprintf(pNum,"%d", intensity); + Tcl_DStringAppendElement(&result,pNum); + Tcl_DStringEndSublist(&result); + } + } + } + Tcl_DStringEndSublist(&result); + SCWrite(pCon,Tcl_DStringValue(&result),eValue); + Tcl_DStringFree(&result); + return 1; + } + else if(strcmp(argv[1],"cog") == 0) /* COG calculation */ + { + if(argc < 5) + { + sprintf(pBueffel,"ERROR: insufficient number of arguments to %s.cog", + argv[0]); + SCWrite(pCon,pBueffel,eError); + return 0; + } + if(!checkHM(&pHM, pSics,pCon, argv[2],iDim)) + { + return 0; + } + if(Tcl_GetInt(pSics->pTcl,argv[3],&i)!= TCL_OK) + { + sprintf(pBueffel,"ERROR: expected number, got %s",argv[2]); + SCWrite(pCon,pBueffel,eError); + return 0; + } + if(Tcl_GetInt(pSics->pTcl,argv[4],&j)!= TCL_OK) + { + sprintf(pBueffel,"ERROR: expected number, got %s",argv[2]); + SCWrite(pCon,pBueffel,eError); + return 0; + } + Tcl_DStringInit(&result); + Tcl_DStringStartSublist(&result); + iData = GetHistogramPointer(pHM,pCon); + window = (int)ObVal(self->pParam,COGWINDOW); + iRet = calculateCOG(iData,iDim[0], iDim[1], &i, &j, &intensity, + window, ObVal(self->pParam,COGCONTOUR)); + if(!iRet) + { + SCWrite(pCon,"ERROR: no intensity in data",eError); + return 0; + } + sprintf(pNum,"%d", i); + Tcl_DStringAppendElement(&result,pNum); + sprintf(pNum,"%d", j); + Tcl_DStringAppendElement(&result,pNum); + sprintf(pNum,"%d", intensity); + Tcl_DStringAppendElement(&result,pNum); + Tcl_DStringEndSublist(&result); + SCWrite(pCon,Tcl_DStringValue(&result),eValue); + Tcl_DStringFree(&result); + return 1; + } + else + { + /* we are handling one of the parameter commands + */ + if(argc > 2) /* set case */ + { + if(Tcl_GetDouble(pSics->pTcl,argv[2],&dVal) != TCL_OK) + { + sprintf(pBueffel,"ERROR: expected number, got %s",argv[2]); + SCWrite(pCon,pBueffel,eError); + return 0; + } + return ObParSet(self->pParam,argv[0],argv[1],(float)dVal,pCon); + } + else /* read case */ + { + ob = ObParFind(self->pParam,argv[1]); + if(!ob) + { + sprintf(pBueffel,"ERROR: parameter %s not found or wrong command", + argv[1]); + SCWrite(pCon,pBueffel,eError); + return 0; + } + sprintf(pBueffel,"%s.%s = %f",argv[0],argv[1],ob->fVal); + SCWrite(pCon,pBueffel,eError); + return 1; + } + } + /* + not reached + */ + assert(0); + return 0; +} +/*-------------------------------------------------------------------*/ +static void KillLoMax(void *pData) +{ + pLoMax self = (pLoMax)pData; + if(!self) + return; + + if(self->pDes) + DeleteDescriptor(self->pDes); + if(self->pParam) + ObParDelete(self->pParam); + free(self); +} +/*-------------------------------------------------------------------*/ +int LoMaxFactory(SConnection *pCon, SicsInterp *pSics, + void *pData, int argc, char *argv[]) +{ + pLoMax pNew = NULL; + + if(argc < 2) + { + SCWrite(pCon,"ERROR: Insufficient number of arguments to LoMaxfactory", + eError); + return 0; + } + + pNew = (pLoMax)malloc(sizeof(LoMax)); + if(!pNew) + { + SCWrite(pCon,"ERROR: out of memory creating local maximum searcher", + eError); + return 0; + } + memset(pNew,0,sizeof(LoMax)); + + /* + create Descriptor + */ + pNew->pDes = CreateDescriptor("Local Maximum Detector"); + if(!pNew->pDes) + { + KillLoMax(pNew); + SCWrite(pCon,"ERROR: out of memory creating local maximum searcher", + eError); + return 0; + } + + /* + create and install parameters + */ + pNew->pParam = ObParCreate(5); + if(!pNew->pParam) + { + KillLoMax(pNew); + SCWrite(pCon,"ERROR: out of memory creating local maximum searcher", + eError); + return 0; + } + ObParInit(pNew->pParam,WINDOW,"window",10,usUser); + ObParInit(pNew->pParam,THRESHOLD,"threshold",30,usUser); + ObParInit(pNew->pParam,STEEPNESS,"steepness",5,usUser); + ObParInit(pNew->pParam,COGWINDOW,"cogwindow",50,usUser); + ObParInit(pNew->pParam,COGCONTOUR,"cogcontour",.2,usUser); + + return AddCommand(pSics, + argv[1], + LoMaxAction, + KillLoMax, + pNew); +} + diff --git a/lomax.h b/lomax.h new file mode 100644 index 00000000..b27c338a --- /dev/null +++ b/lomax.h @@ -0,0 +1,37 @@ + +/*------------------------------------------------------------------------- + L o c a l M a x i m u m + + This is a module for searching a local maximum in a 2D histogram as + collected at TRICS. + + copyright: see copyright.h + + Mark Koennecke, November 2001 +-------------------------------------------------------------------------*/ +#ifndef LOCALMAXIMUM +#define LOCALMAXIMUM +#include "obpar.h" + + typedef struct __LOMAX *pLoMax; + + + int LoMaxAction(SConnection *pCon, SicsInterp *pSics, + void *pData, int argc, char *argv[]); + + int LoMaxFactory(SConnection *pCon, SicsInterp *pSics, + void *pData, int argc, char *argv[]); + + int testLocalMaximum(int *iData, int xsize, int ysize, + int i, int j, + int window, int threshold, int steepness, + int *intensity); + int calculateCOG(int *iData, int xsize, int ysize, + int *i, int *j, int *intensity, + int cogWindow, + float contour); + + + + +#endif diff --git a/lomax.i b/lomax.i new file mode 100644 index 00000000..fd803f79 --- /dev/null +++ b/lomax.i @@ -0,0 +1,19 @@ + +/*------------------------------------------------------------------------- + L o c a l M a x i m u m + + This is a module for searching a local maximum in a 2D histogram as + collected at TRICS. These are internal definitions only and are + supposed to be included only into lomax.c + + copyright: see copyright.h + + Mark Koennecke, November 2001 +-------------------------------------------------------------------------*/ + + typedef struct __LOMAX { + pObjectDescriptor pDes; + ObPar *pParam; + }LoMax; + + diff --git a/lomax.w b/lomax.w new file mode 100644 index 00000000..d6684646 --- /dev/null +++ b/lomax.w @@ -0,0 +1,104 @@ +\subsection{Local Maximum Search} +This module implements searching for local maxima in a 2D histogram of +data. It is mainly intended for use at TRICS where initially +reflections need to be searched for the determination of a UB +matrix. The local maximum algorithm is a little more sophisticated +then commonly used in that the window of the search is user definable, +and a valid maximum has to have a certain steepness against the +bordering points in order to be valid. Another validity test is a +threshold which must be reached. All these enhancements should reduce +the number of false maxima which may turn up as the result of horrible +measurement statistics. + +The position and intensity of a local maximum found can be refined further + with a center of gravity calculation within a given window including all + data points above a certain threshold expressed as a fraction of the local + maximums intensity. + +Again this module has a little data structure: + +@d lomaxdata @{ + typedef struct __LOMAX { + pObjectDescriptor pDes; + ObPar *pParam; + }LoMax; +@} +The fields are: +\begin{description} +\item[pDes] The standard object descriptor. +\item[pParam] The parameters of this object. +\end{description} + +The exported interface is simple: just the necessary functions for +representing the object in the interpreter and a function which allows +to use the local maxima search from another module. If the latter is +ever needed. + +@d lomaxint @{ + int LoMaxAction(SConnection *pCon, SicsInterp *pSics, + void *pData, int argc, char *argv[]); + + int LoMaxFactory(SConnection *pCon, SicsInterp *pSics, + void *pData, int argc, char *argv[]); + + int testLocalMaximum(int *iData, int xsize, int ysize, + int i, int j, + int window, int threshold, int steepness, + int *intensity); + int calculateCOG(int *iData, int xsize, int ysize, + int *i, int *j, int *intensity, + int cogWindow, + float contour); + +@} + +testLocalMaxima checks if point i, j is a valid local maximum. It +returns true (1) if this is the case, or false (0) if this is not the +case. The window parameter describes how much points shall be searched +in each direction. + +calculateCOG calculates the center of gravity for point i, j,. Points within a + sqaure window of size cogWindow are considered. Only points whose + intensity obeys I \gt I(i,j)*countour are used for the COG +calculation. The position of the maximum (i,j) and its intensity are +updated by calculateCOG. + + +@o lomax.h @{ +/*------------------------------------------------------------------------- + L o c a l M a x i m u m + + This is a module for searching a local maximum in a 2D histogram as + collected at TRICS. + + copyright: see copyright.h + + Mark Koennecke, November 2001 +-------------------------------------------------------------------------*/ +#ifndef LOCALMAXIMUM +#define LOCALMAXIMUM +#include "obpar.h" + + typedef struct __LOMAX *pLoMax; + +@ + + +#endif +@} + +@o lomax.i @{ +/*------------------------------------------------------------------------- + L o c a l M a x i m u m + + This is a module for searching a local maximum in a 2D histogram as + collected at TRICS. These are internal definitions only and are + supposed to be included only into lomax.c + + copyright: see copyright.h + + Mark Koennecke, November 2001 +-------------------------------------------------------------------------*/ +@ + +@} diff --git a/nextrics.c b/nextrics.c index 43bcf286..63143d29 100644 --- a/nextrics.c +++ b/nextrics.c @@ -395,8 +395,8 @@ name of hkl object holding crystallographic information NXputdata(hfil,lData); NXsetdimname(hfil,0,"frame_x"); NXsetdimname(hfil,1,"frame_y"); - */ NXclosedata(hfil); + */ fVal = fTTheta + self->hm1off; NXDputalias(hfil,self->pDict,"frame2theta",&fVal); @@ -424,8 +424,8 @@ name of hkl object holding crystallographic information NXputdata(hfil,lData); NXsetdimname(hfil,0,"frame_x"); NXsetdimname(hfil,1,"frame_y"); - */ NXclosedata(hfil); + */ fVal = fTTheta + self->hm2Off; NXDputalias(hfil,self->pDict,"frame2theta",&fVal); @@ -454,9 +454,9 @@ name of hkl object holding crystallographic information NXputdata(hfil,lData); NXsetdimname(hfil,0,"frame_x"); NXsetdimname(hfil,1,"frame_y"); + NXclosedata(hfil); */ - NXclosedata(hfil); fVal = fTTheta + self->hm3Off; NXDputalias(hfil,self->pDict,"frame2theta",&fVal); diff --git a/ofac.c b/ofac.c index 968b47df..a60e87b9 100644 --- a/ofac.c +++ b/ofac.c @@ -105,6 +105,7 @@ #include "swmotor.h" #include "hmcontrol.h" #include "rs232controller.h" +#include "lomax.h" /*----------------------- Server options creation -------------------------*/ static int IFServerOption(SConnection *pCon, SicsInterp *pSics, void *pData, int argc, char *argv[]) @@ -285,6 +286,7 @@ AddCommand(pInter,"MakeSWMotor",MakeSWMotor,NULL,NULL); AddCommand(pInter,"MakeHMControl",MakeHMControl,NULL,NULL); AddCommand(pInter,"MakeRS232Controller",RS232Factory,NULL,NULL); + AddCommand(pInter,"MakeMaxDetector",LoMaxFactory,NULL,NULL); } /*---------------------------------------------------------------------------*/ static void KillIniCommands(SicsInterp *pSics) @@ -342,6 +344,7 @@ RemoveCommand(pSics,"MakeSWMotor"); RemoveCommand(pSics,"MakeHMControl"); RemoveCommand(pSics,"MakeRS232Controller"); + RemoveCommand(pSics,"MakeMaxDetector"); } diff --git a/peaksearch.tcl b/peaksearch.tcl new file mode 100644 index 00000000..c5b6423b --- /dev/null +++ b/peaksearch.tcl @@ -0,0 +1,208 @@ +#--------------------------------------------------------------------------- +# peaksearch a peak search utility for TRICS using the PSD detectors. +# +# Mark Koennecke, November 2001 +#--------------------------------------------------------------------------- + +proc initPeakSearch {} { + VarMake ps.phiStart Float User + ps.phiStart 0 + VarMake ps.phiEnd Float User + ps.phiEnd 180 + VarMake ps.phiStep Float User + ps.phiStep 3. + + VarMake ps.chiStart Float User + ps.chiStart 0 + VarMake ps.chiEnd Float User + ps.chiEnd 180 + VarMake ps.chiStep Float User + ps.chiStep 12. + + VarMake ps.omStart Float User + ps.omStart 0 + VarMake ps.omEnd Float User + ps.omEnd 30 + VarMake ps.omStep Float User + ps.omStep 3. + + VarMake ps.sttStart Float User + ps.sttStart 5 + VarMake ps.sttEnd Float User + ps.sttEnd 70 + VarMake ps.sttStep Float User + ps.sttStep 3. + + VarMake ps.threshold Int User + ps.threshold 30 + VarMake ps.steepness Int User + ps.steepness 3 + VarMake ps.window Int User + ps.window 7 + VarMake ps.cogwindow Int User + ps.cogwindow 60 + VarMake ps.cogcontour Float User + ps.cogcontour .2 + + VarMake ps.countmode Text User + ps.countmode monitor + + VarMake ps.preset Float User + ps.preset 1000 + + Publish ps.phirange User + Publish ps.chirange User + Publish ps.omrange User + Publish ps.sttrange User + Publish ps.countpar User + Publish ps.maxpar User +#------- these are for debugging only! + Publish checkomega User + Publish optimizepeak User +} +#-------------------------------------------------------------------------- +proc ps.phirange args { + if { [llength $args] >= 3 } { + ps.phiStart [lindex $args 0] + ps.phiEnd [lindex $args 1] + ps.phiStep [lindex $args 2] + } + clientput "Peak Search Phi Range:" + return [format " Start = %6.2f, End = %6.2f, Step = %6.2f" \ + [SplitReply [ps.phiStart]] [SplitReply [ps.phiEnd]] \ + [SplitReply [ps.phiStep]]] +} +#-------------------------------------------------------------------------- +proc ps.chirange args { + if { [llength $args] >= 3 } { + ps.chiStart [lindex $args 0] + ps.chiEnd [lindex $args 1] + ps.chiStep [lindex $args 2] + } + clientput "Peak Search Chi Range:" + return [format " Start = %6.2f, End = %6.2f, Step = %6.2f" \ + [SplitReply [ps.chiStart]] [SplitReply [ps.chiEnd]] \ + [SplitReply [ps.chiStep]]] +} +#-------------------------------------------------------------------------- +proc ps.omrange args { + if { [llength $args] >= 3 } { + ps.omStart [lindex $args 0] + ps.omEnd [lindex $args 1] + ps.omStep [lindex $args 2] + } + clientput "Peak Search Omega Range:" + return [format " Start = %6.2f, End = %6.2f, Step = %6.2f" \ + [SplitReply [ps.omStart]] [SplitReply [ps.omEnd]] \ + [SplitReply [ps.omStep]]] +} +#-------------------------------------------------------------------------- +proc ps.sttrange args { + if { [llength $args] >= 3 } { + ps.sttStart [lindex $args 0] + ps.sttEnd [lindex $args 1] + ps.sttStep [lindex $args 2] + } + clientput "Peak Search Two Theta Range:" + return [format " Start = %6.2f, End = %6.2f, Step = %6.2f" \ + [SplitReply [ps.sttStart]] [SplitReply [ps.sttEnd]] \ + [SplitReply [ps.sttStep]]] +} +#------------------------------------------------------------------------- +proc ps.countpar args { + if { [llength $args] >= 2 } { + if { [catch {counter setmode [lindex $args 0]} msg] != 0} { + error "ERROR: Invalid countmode specified" + } + ps.countmode [lindex $args 0] + ps.preset [lindex $args 1] + } + clientput "Peak Search Count Parameters:" + return [format " Mode = %s, Preset = %12.2f" \ + [SplitReply [ps.countmode]] [SplitReply [ps.preset]]] +} +#------------------------------------------------------------------------- +proc ps.maxpar args { + if { [llength $args] >= 5 } { + ps.window [lindex $args 0] + ps.threshold [lindex $args 1] + ps.steepness [lindex $args 2] + ps.cogwindow [lindex $args 3] + ps.cogcontour [lindex $args 4] + } + clientput "Peak Search Maximum Detection Parameters:" + set t1 [format " Window = %d, Threshold = %d, Steepness = %d" \ + [SplitReply [ps.window]] [SplitReply [ps.threshold]] \ + [SplitReply [ps.steepness] ]] + set t2 [format " COGWindow = %d, COGcontour = %f " \ + [SplitReply [ps.cogwindow]] [SplitReply [ps.cogcontour]]] + return [format "%s\n%s" $t1 $t2] +} +#------------------------------------------------------------------------ +proc checknewomega {omega maxIntensity} { + if {[catch {drive $omega} msg] != 0} { + error $msg + } + if {[catch {hmc start [SplitReply [ps.preset]] [SplitReply \ + [ps.countmode]]} msg] != 0} { + error $msg + } + if {[catch {Success} msg] != 0} { + error $msg + } + if { [catch {lomax cog $hm $x $y} result] != 0} { + error "Failed to calculate COG" + } + if {[lindex $result 2] > $maxIntensity } { + return $result + } else { + return 0 + } +} +#------------------------------------------------------------------------ +proc optimizepeak {hm x y} { + if { [catch {lomax cog $hm $x $y} result] != 0} { + error "Failed to calculate COG" + } + set xMax [lindex $result 0] + set ymax [lindex $result 1] + set maxIntensity [lindex $result 2] + set maxOmega [SplitReply [om]] + set startOmega $maxOmega +#--------- move to positive omega until maximum found + while {1} { + set newOm [expr [SplitReply [om]] + [SplitReply [ps.omStep]]] + if {[catch {checknewomega $newOm $maxIntensity} result] != 0} { + error $result + } + if {$result != 0} { + set xMax [lindex $result 0] + set yMax [lindex $result 1] + set maxIntensity [lindex $result 2] + set maxOmega [SplitReply [om]] + } else { + break + } + } +#--------- if maxOmega is still startOmega then we were on the right +# side of the peak. In this case try to find the maximum in +# negative direction + if {$maxOmega == $startOmega} { + while {1} { + set newOm [expr [SplitReply [om]] - [SplitReply [ps.omStep]]] + if {[catch {checknewomega $newOm $maxIntensity} result] != 0} { + error $result + } + if {$result != 0} { + set xMax [lindex $result 0] + set yMax [lindex $result 1] + set maxIntensity [lindex $result 2] + set maxOmega [SplitReply [om]] + } else { + break + } + } + } +#----------- print the results we have found +} + diff --git a/sicsstatus.tcl b/sicsstatus.tcl index 1e7e03bd..8d6f7ef1 100644 --- a/sicsstatus.tcl +++ b/sicsstatus.tcl @@ -1,30 +1,47 @@ -a5l.length 80.000000 -flightpathlength 0.000000 -flightpathlength setAccess 1 -flightpath 0.000000 -flightpath setAccess 1 -delay 2500.000000 -delay setAccess 1 -hm CountMode timer -hm preset 100.000000 -hm genbin 120.000000 35.000000 512 -hm init -datafile focus-1001848.hdf -datafile setAccess 3 +ps.preset 20.000000 +ps.preset setAccess 2 +ps.countmode timer +ps.countmode setAccess 2 +ps.cogcontour 0.200000 +ps.cogcontour setAccess 2 +ps.cogwindow 60 +ps.cogwindow setAccess 2 +ps.window 7 +ps.window setAccess 2 +ps.steepness 2 +ps.steepness setAccess 2 +ps.threshold 50 +ps.threshold setAccess 2 +ps.sttstep 5.000000 +ps.sttstep setAccess 2 +ps.sttend 50.000000 +ps.sttend setAccess 2 +ps.sttstart 10.000000 +ps.sttstart setAccess 2 +ps.omstep 5.000000 +ps.omstep setAccess 2 +ps.omend 20.000000 +ps.omend setAccess 2 +ps.omstart 10.000000 +ps.omstart setAccess 2 +ps.chistep 10.000000 +ps.chistep setAccess 2 +ps.chiend 150.000000 +ps.chiend setAccess 2 +ps.chistart 10.000000 +ps.chistart setAccess 2 +ps.phistep 5.000000 +ps.phistep setAccess 2 +ps.phiend 120.000000 +ps.phiend setAccess 2 +ps.phistart 10.000000 +ps.phistart setAccess 2 +hm3 CountMode timer +hm3 preset 10.000000 hm2 CountMode timer hm2 preset 5.000000 hm1 CountMode timer hm1 preset 5.000000 -dbfile UNKNOWN -dbfile setAccess 2 -# Motor th -th SoftZero 0.000000 -th SoftLowerLim -120.000000 -th SoftUpperLim 120.000000 -th Fixed -1.000000 -th sign 1.000000 -th InterruptMode 0.000000 -th AccessCode 2.000000 #Crystallographic Settings hkl lambda 1.179000 hkl setub 0.016169 0.011969 0.063195 -0.000545 0.083377 -0.009117 -0.162051 0.000945 0.006312 @@ -148,8 +165,6 @@ twotheta InterruptMode 0.000000 twotheta AccessCode 2.000000 lastscancommand cscan a4 10. .1 10 5 lastscancommand setAccess 2 -banana CountMode timer -banana preset 100.000000 sample_mur 0.000000 sample_mur setAccess 2 email UNKNOWN