- Timestamp:
- Sep 13, 2019 3:19:52 PM (3 years ago)
- Location:
- sans/Dev/branches/elliptical_averaging/NCNR_User_Procedures/Reduction/SANS
- Files:
-
- 2 edited
-
AvgGraphics.ipf (modified) (4 diffs)
-
RectAnnulAvg.ipf (modified) (2 diffs)
Legend:
- Unmodified
- Added
- Removed
-
sans/Dev/branches/elliptical_averaging/NCNR_User_Procedures/Reduction/SANS/AvgGraphics.ipf
r1207 r1208 327 327 GroupBox ann,pos={148,44},size={143,84},title="Annular" 328 328 GroupBox rect,pos={7,133},size={134,71},title="Rectangular" 329 GroupBox sect,pos={148,133},size={143,71},title="Sector "329 GroupBox sect,pos={148,133},size={143,71},title="Sector/Elliptical" 330 330 GroupBox sect_rect,pos={7,44},size={134,84},title="Sector/Rectangular" 331 331 PopupMenu av_choice,pos={61,7},size={144,20},proc=AvTypePopMenuProc,title="AverageType" … … 476 476 Return 0 //exit the Draw routine 477 477 Endif 478 479 // Elliptical 480 if(cmpstr(av_type,"Elliptical")==0) 481 // 482 // TODO: Write the elliptical averaging design 483 // 484 485 //do the elliptical drawing 486 sdd=reals[18] 487 lam=reals[26] 488 489 sdd *=100 //convert to cm 490 //find the radius (from the y-direction) 491 variable RAxes = NumberByKey("RATIOAXES",drawStr,"=",";") 492 493 Make/O/N=(128,128) ellipsewave 494 variable b = 64 495 variable a = RAxes * b 496 DrawAnEllipsoid("ellipsewave",x0,y0,a,b,phi) 497 498 Return 0 //exit the Draw routine 499 EndIf 478 500 479 501 //else sector or rectangular - draw the lines … … 795 817 End 796 818 819 820 // This creates an ellipsoid 821 // INPUTS: 822 // wname - name of xywave 823 // (xo,yo) - center position of ellipse 824 // (a, b) - elliptical dimesions to above 825 // phi - angle of rotation (counterclockwise) 826 // OUTPUTS: 827 // make all changes directly in xywave 828 Function DrawAnEllipsoid(wname,xo,yo,a,b,phi) 829 string wname 830 variable xo,yo,a,b,phi 831 int ii 832 833 wave xywave = $wname 834 835 make/FREE/n=361 urx,ury 836 make/n=361 xw,yw = 0 837 838 variable cosalpha = cos(phi*Pi/180), sinalpha = sin(phi*Pi/180) 839 840 urx = a*(cos(p*Pi/180)) 841 ury = b*(sin(p*Pi/180)) 842 843 xw = xo + cosalpha* urx - sinalpha*ury 844 yw = yo + sinalpha*urx + cosalpha*ury 845 846 // FIXME: No single point mapping in Igor - need to do this automatically 847 848 xywave[xw,yw] = 1 849 850 rr=0 851 gg=50000 852 bb=0 853 854 //SetDrawEnv/W=SANS_Data xcoord= bottom,ycoord= left,linefgc= (rr,gg,bb),linethick= (thick),fillpat=0 855 //AppendImage/W=SANS_Data ellipsewave 856 857 Return (0) 858 End 859 797 860 //reads the value from the panel and updates the global string 798 861 // redraws the appropriate line on the image … … 866 929 867 930 String/G root:myGlobals:Drawing:gDrawInfoStr = newStr 931 932 //redraw the angles 933 MasterAngleDraw() 934 End 935 936 //reads the value from the panel and updates the global string 937 // redraws the appropriate line on the image 938 Function AxisRatioSetVarProc(ctrlName,varNum,varStr,varName) : SetVariableControl 939 String ctrlName 940 Variable varNum 941 String varStr 942 String varName 943 944 ControlInfo/W=Average_Panel RAxes_p 945 Variable val = V_Value 946 SVAR tempStr = root:myGlobals:Drawing:gDrawRAxes 947 String newStr = ReplaceNumberByKey("RATIOAXES", tempStr, val, "=", ";") 948 949 String/G root:myGlobals:Drawing:gDrawRAxes = newStr 868 950 869 951 //redraw the angles -
sans/Dev/branches/elliptical_averaging/NCNR_User_Procedures/Reduction/SANS/RectAnnulAvg.ipf
r1207 r1208 706 706 lambda = reals[26] 707 707 708 Variable qc = NumberByKey("QCENTER",keyListStr,"=",";")709 Variable nw = NumberByKey("QDELTA",keyListStr,"=",";")708 Variable phi = NumberByKey("PHI",keyListStr,"=",";") 709 Variable ratio = NumberByKey("RATIOAXES",keyListStr,"=",";") 710 710 711 711 dr = 1 //minimum annulus width, keep this fixed at one … … 713 713 nphi = numPhiSteps //number of anular sectors is set by users 714 714 715 rc = 2*dtdist*asin(qc*lambda/4/Pi) //in mm716 delr = nw*sx/2717 rlo = rc-delr718 rhi = rc + delr719 dphi = 360/nphi720 721 /// data wave is data in the current folder which was set at the top of the function722 Wave data=$(destPath + ":data")723 //Check for the existence of the mask, if not, make one (local to this folder) that is null724 725 if(WaveExists($"root:Packages:NIST:MSK:data") == 0)726 Print "There is no mask file loaded (WaveExists)- the data is not masked"727 Make/O/N=(pixelsX,pixelsY) $(destPath + ":mask")728 WAVE mask = $(destPath + ":mask")729 mask = 0730 else731 Wave mask=$"root:Packages:NIST:MSK:data"732 Endif733 734 rcentr = 150 //pixels within rcentr of beam center are broken into 9 parts735 // values for error if unable to estimate value736 //large_num = 1e10737 large_num = 1 //1e10 value (typically sig of last data point) plots poorly, arb set to 1738 small_num = 1e-10739 740 // output wave are expected to exist (?) initialized to zero, what length?741 // 300 points on VAX ---742 Variable wavePts=500743 Make/O/N=(wavePts) $(destPath + ":phival"),$(destPath + ":aveint")744 Make/O/N=(wavePts) $(destPath + ":ncells"),$(destPath + ":sig"),$(destPath + ":sigave")745 WAVE phival = $(destPath + ":phival")746 WAVE aveint = $(destPath + ":aveint")747 WAVE ncells = $(destPath + ":ncells")748 WAVE sig = $(destPath + ":sig")749 WAVE sigave = $(destPath + ":sigave")750 751 phival = 0752 aveint = 0753 ncells = 0754 sig = 0755 sigave = 0756 757 dtdis2 = dtdist^2758 nq = 1759 xoffst=0760 //distance of beam center from detector center761 xbm = FX(x0,sx3,xcenter,sx)762 ybm = FY(y0,sy3,ycenter,sy)763 764 //BEGIN AVERAGE **********765 Variable xi,xd,x,y,yd,yj,nd,fd,nd2,iphi,ntotal,var766 Variable jj,data_pixel,xx,yy,ll,kk,rij,phiij,avesq,aveisq767 768 // IGOR arrays are indexed from [0][0], FORTAN from (1,1) (and the detector too)769 // loop index corresponds to FORTRAN (old code)770 // and the IGOR array indices must be adjusted (-1) to the correct address771 ntotal = 0772 ii=1773 do774 xi = ii775 xd = FX(xi,sx3,xcenter,sx)776 x = xoffst + xd -xbm //x and y are in mm777 778 jj = 1779 do780 data_pixel = data[ii-1][jj-1] //assign to local variable781 yj = jj782 yd = FY(yj,sy3,ycenter,sy)783 y = yd - ybm784 if(!(mask[ii-1][jj-1])) //masked pixels = 1, skip if masked (this way works...)785 nd = 1786 fd = 1787 if( (abs(x) > rcentr) || (abs(y) > rcentr)) //break pixel into 9 equal parts788 nd = 3789 fd = 2790 Endif791 nd2 = nd^2792 ll = 1 //"el-el" loop index793 do794 xx = x + (ll - fd)*sx/3795 kk = 1796 do797 yy = y + (kk - fd)*sy/3798 //test to see if center of pixel (i,j) lies in annulus799 rij = sqrt(x*x + y*y)/dr + 1.001800 //check whether pixel lies within width band801 if((rij > rlo) && (rij < rhi))802 //in the annulus, do something803 if (yy >= 0)804 //phiij is in degrees805 phiij = atan2(yy,xx)*180/Pi //0 to 180 deg806 else807 phiij = 360 + atan2(yy,xx)*180/Pi //180 to 360 deg808 Endif809 if (phiij > (360-0.5*dphi))810 phiij -= 360811 Endif812 iphi = trunc(phiij/dphi + 1.501)813 aveint[iphi-1] += 9*data_pixel/nd2814 sig[iphi-1] += 9*data_pixel*data_pixel/nd2815 ncells[iphi-1] += 9/nd2816 ntotal += 9/nd2817 Endif //check if in annulus818 kk+=1819 while(kk<=nd)820 ll += 1821 while(ll<=nd)822 Endif //masked pixel check823 jj += 1824 while (jj<=pixelsY)825 ii += 1826 while(ii<=pixelsX) //end of the averaging827 828 //compute phi-values and errors829 830 ntotal /=9831 832 kk = 1833 do834 phival[kk-1] = dphi*(kk-1)835 if(ncells[kk-1] != 0)836 aveint[kk-1] = aveint[kk-1]/ncells[kk-1]837 avesq = aveint[kk-1]*aveint[kk-1]838 aveisq = sig[kk-1]/ncells[kk-1]839 var = aveisq - avesq840 if (var <=0 )841 sig[kk-1] = 0842 sigave[kk-1] = 0843 ncells[kk-1] /=9844 else845 if(ncells[kk-1] > 9)846 sigave[kk-1] = sqrt(9*var/(ncells[kk-1]-9))847 sig[kk-1] = sqrt( abs(aveint[kk-1])/(ncells[kk-1]/9) )848 ncells[kk-1] /=9849 else850 sig[kk-1] = 0851 sigave[kk-1] = 0852 ncells[kk-1] /=9853 Endif854 Endif855 Endif856 kk+=1857 while(kk<=nphi)858 859 // data waves were defined as 200 points (=wavePts), but now have less than that (nphi) points860 // use DeletePoints to remove junk from end of waves861 Variable startElement,numElements862 startElement = nphi863 numElements = wavePts - startElement864 DeletePoints startElement,numElements, phival,aveint,ncells,sig,sigave865 866 //////////////end of VAX Phibin.for867 868 //angle dependent transmission correction is not done in phiave869 Ann_1D_Graph(aveint,phival,sigave)715 // rc = 2*dtdist*asin(qc*lambda/4/Pi) //in mm 716 // delr = nw*sx/2 717 // rlo = rc-delr 718 // rhi = rc + delr 719 // dphi = 360/nphi 720 // 721 // /// data wave is data in the current folder which was set at the top of the function 722 // Wave data=$(destPath + ":data") 723 // //Check for the existence of the mask, if not, make one (local to this folder) that is null 724 // 725 // if(WaveExists($"root:Packages:NIST:MSK:data") == 0) 726 // Print "There is no mask file loaded (WaveExists)- the data is not masked" 727 // Make/O/N=(pixelsX,pixelsY) $(destPath + ":mask") 728 // WAVE mask = $(destPath + ":mask") 729 // mask = 0 730 // else 731 // Wave mask=$"root:Packages:NIST:MSK:data" 732 // Endif 733 // 734 // rcentr = 150 //pixels within rcentr of beam center are broken into 9 parts 735 // // values for error if unable to estimate value 736 // //large_num = 1e10 737 // large_num = 1 //1e10 value (typically sig of last data point) plots poorly, arb set to 1 738 // small_num = 1e-10 739 // 740 // // output wave are expected to exist (?) initialized to zero, what length? 741 // // 300 points on VAX --- 742 // Variable wavePts=500 743 // Make/O/N=(wavePts) $(destPath + ":phival"),$(destPath + ":aveint") 744 // Make/O/N=(wavePts) $(destPath + ":ncells"),$(destPath + ":sig"),$(destPath + ":sigave") 745 // WAVE phival = $(destPath + ":phival") 746 // WAVE aveint = $(destPath + ":aveint") 747 // WAVE ncells = $(destPath + ":ncells") 748 // WAVE sig = $(destPath + ":sig") 749 // WAVE sigave = $(destPath + ":sigave") 750 // 751 // phival = 0 752 // aveint = 0 753 // ncells = 0 754 // sig = 0 755 // sigave = 0 756 // 757 // dtdis2 = dtdist^2 758 // nq = 1 759 // xoffst=0 760 // //distance of beam center from detector center 761 // xbm = FX(x0,sx3,xcenter,sx) 762 // ybm = FY(y0,sy3,ycenter,sy) 763 // 764 // //BEGIN AVERAGE ********** 765 // Variable xi,xd,x,y,yd,yj,nd,fd,nd2,iphi,ntotal,var 766 // Variable jj,data_pixel,xx,yy,ll,kk,rij,phiij,avesq,aveisq 767 // 768 // // IGOR arrays are indexed from [0][0], FORTAN from (1,1) (and the detector too) 769 // // loop index corresponds to FORTRAN (old code) 770 // // and the IGOR array indices must be adjusted (-1) to the correct address 771 // ntotal = 0 772 // ii=1 773 // do 774 // xi = ii 775 // xd = FX(xi,sx3,xcenter,sx) 776 // x = xoffst + xd -xbm //x and y are in mm 777 // 778 // jj = 1 779 // do 780 // data_pixel = data[ii-1][jj-1] //assign to local variable 781 // yj = jj 782 // yd = FY(yj,sy3,ycenter,sy) 783 // y = yd - ybm 784 // if(!(mask[ii-1][jj-1])) //masked pixels = 1, skip if masked (this way works...) 785 // nd = 1 786 // fd = 1 787 // if( (abs(x) > rcentr) || (abs(y) > rcentr)) //break pixel into 9 equal parts 788 // nd = 3 789 // fd = 2 790 // Endif 791 // nd2 = nd^2 792 // ll = 1 //"el-el" loop index 793 // do 794 // xx = x + (ll - fd)*sx/3 795 // kk = 1 796 // do 797 // yy = y + (kk - fd)*sy/3 798 // //test to see if center of pixel (i,j) lies in annulus 799 // rij = sqrt(x*x + y*y)/dr + 1.001 800 // //check whether pixel lies within width band 801 // if((rij > rlo) && (rij < rhi)) 802 // //in the annulus, do something 803 // if (yy >= 0) 804 // //phiij is in degrees 805 // phiij = atan2(yy,xx)*180/Pi //0 to 180 deg 806 // else 807 // phiij = 360 + atan2(yy,xx)*180/Pi //180 to 360 deg 808 // Endif 809 // if (phiij > (360-0.5*dphi)) 810 // phiij -= 360 811 // Endif 812 // iphi = trunc(phiij/dphi + 1.501) 813 // aveint[iphi-1] += 9*data_pixel/nd2 814 // sig[iphi-1] += 9*data_pixel*data_pixel/nd2 815 // ncells[iphi-1] += 9/nd2 816 // ntotal += 9/nd2 817 // Endif //check if in annulus 818 // kk+=1 819 // while(kk<=nd) 820 // ll += 1 821 // while(ll<=nd) 822 // Endif //masked pixel check 823 // jj += 1 824 // while (jj<=pixelsY) 825 // ii += 1 826 // while(ii<=pixelsX) //end of the averaging 827 // 828 // //compute phi-values and errors 829 // 830 // ntotal /=9 831 // 832 // kk = 1 833 // do 834 // phival[kk-1] = dphi*(kk-1) 835 // if(ncells[kk-1] != 0) 836 // aveint[kk-1] = aveint[kk-1]/ncells[kk-1] 837 // avesq = aveint[kk-1]*aveint[kk-1] 838 // aveisq = sig[kk-1]/ncells[kk-1] 839 // var = aveisq - avesq 840 // if (var <=0 ) 841 // sig[kk-1] = 0 842 // sigave[kk-1] = 0 843 // ncells[kk-1] /=9 844 // else 845 // if(ncells[kk-1] > 9) 846 // sigave[kk-1] = sqrt(9*var/(ncells[kk-1]-9)) 847 // sig[kk-1] = sqrt( abs(aveint[kk-1])/(ncells[kk-1]/9) ) 848 // ncells[kk-1] /=9 849 // else 850 // sig[kk-1] = 0 851 // sigave[kk-1] = 0 852 // ncells[kk-1] /=9 853 // Endif 854 // Endif 855 // Endif 856 // kk+=1 857 // while(kk<=nphi) 858 // 859 // // data waves were defined as 200 points (=wavePts), but now have less than that (nphi) points 860 // // use DeletePoints to remove junk from end of waves 861 // Variable startElement,numElements 862 // startElement = nphi 863 // numElements = wavePts - startElement 864 // DeletePoints startElement,numElements, phival,aveint,ncells,sig,sigave 865 // 866 // //////////////end of VAX Phibin.for 867 // 868 // //angle dependent transmission correction is not done in phiave 869 // Ann_1D_Graph(aveint,phival,sigave) 870 870 871 871 //get rid of the default mask, if one was created (it is in the current folder)
Note: See TracChangeset
for help on using the changeset viewer.
