 Timestamp:
 Jun 16, 2010 2:10:47 PM (13 years ago)
 File:

 1 edited
Legend:
 Unmodified
 Added
 Removed

sans/Dev/trunk/NCNR_User_Procedures/Common/GaussUtils_v40.ipf
r570 r708 62 62 63 63 64 // tentative pass at 2D resolution smearing 64 //// tentative pass at 2D resolution smearing 65 //// 66 //Structure ResSmear_2D_AAOStruct 67 // Wave coefW 68 // Wave zw //answer 69 // Wave qy // qvalue 70 // Wave qx 71 // Wave qz 72 // Wave sQpl //resolution parallel to Q 73 // Wave sQpp //resolution perpendicular to Q 74 // Wave fs 75 // String info 76 //EndStructure 77 78 // reformat the structure ?? WM Fit compatible 79 //  2D resolution smearing 65 80 // 66 81 Structure ResSmear_2D_AAOStruct 67 82 Wave coefW 68 83 Wave zw //answer 69 Wave qy // qvalue 70 Wave qx 84 Wave xw[2] // qxvalue is [0], qy is xw[1] 71 85 Wave qz 72 Wave s igQx //resolution73 Wave s igQy86 Wave sQpl //resolution parallel to Q 87 Wave sQpp //resolution perpendicular to Q 74 88 Wave fs 75 89 String info … … 1077 1091 1078 1092 // prototype function for 2D smearing routine 1079 Function SANS_2D_ModelAAO_proto(w,zw,xw,yw)1093 ThreadSafe Function SANS_2D_ModelAAO_proto(w,zw,xw,yw) 1080 1094 Wave w,zw,xw,yw 1081 1095 … … 1346 1360 return(0) 1347 1361 end 1362 1363 1364 //// 1365 //// moved from RawWindowHook to here  where the Q calculations are available to 1366 // reduction and analysis 1367 // 1368 1369 //phi is defined from +x axis, proceeding CCW around [0,2Pi] 1370 Threadsafe Function FindPhi(vx,vy) 1371 variable vx,vy 1372 1373 variable phi 1374 1375 phi = atan(vy/vx) //returns a value from pi/2 to pi/2 1376 1377 // special cases 1378 if(vx==0 && vy > 0) 1379 return(pi/2) 1380 endif 1381 if(vx==0 && vy < 0) 1382 return(3*pi/2) 1383 endif 1384 if(vx >= 0 && vy == 0) 1385 return(0) 1386 endif 1387 if(vx < 0 && vy == 0) 1388 return(pi) 1389 endif 1390 1391 1392 if(vx > 0 && vy > 0) 1393 return(phi) 1394 endif 1395 if(vx < 0 && vy > 0) 1396 return(phi + pi) 1397 endif 1398 if(vx < 0 && vy < 0) 1399 return(phi + pi) 1400 endif 1401 if( vx > 0 && vy < 0) 1402 return(phi + 2*pi) 1403 endif 1404 1405 return(phi) 1406 end 1407 1408 1409 //function to calculate the overall qvalue, given all of the necesary trig inputs 1410 //NOTE: detector locations passed in are pixels = 0.5cm real space on the detector 1411 //and are in detector coordinates (1,128) rather than axis values 1412 //the pixel locations need not be integers, reals are ok inputs 1413 //sdd is in meters 1414 //wavelength is in Angstroms 1415 // 1416 //returned magnitude of Q is in 1/Angstroms 1417 // 1418 Function CalcQval(xaxval,yaxval,xctr,yctr,sdd,lam,pixSize) 1419 Variable xaxval,yaxval,xctr,yctr,sdd,lam,pixSize 1420 1421 Variable dx,dy,qval,two_theta,dist 1422 1423 Variable pixSizeX=pixSize 1424 Variable pixSizeY=pixSize 1425 1426 sdd *=100 //convert to cm 1427 dx = (xaxval  xctr)*pixSizeX //delta x in cm 1428 dy = (yaxval  yctr)*pixSizeY //delta y in cm 1429 dist = sqrt(dx^2 + dy^2) 1430 1431 two_theta = atan(dist/sdd) 1432 1433 qval = 4*Pi/lam*sin(two_theta/2) 1434 1435 return qval 1436 End 1437 1438 //calculates just the qvalue in the xdirection on the detector 1439 //input/output is the same as CalcQval() 1440 //ALL inputs are in detector coordinates 1441 // 1442 //NOTE: detector locations passed in are pixel = 0.5cm real space on the Ordela detector 1443 //sdd is in meters 1444 //wavelength is in Angstroms 1445 // 1446 // repaired incorrect qx and qy calculation 3 dec 08 SRK (Lionel and C. Dewhurst) 1447 // now properly accounts for qz 1448 // 1449 Function CalcQX(xaxval,yaxval,xctr,yctr,sdd,lam,pixSize) 1450 Variable xaxval,yaxval,xctr,yctr,sdd,lam,pixSize 1451 1452 Variable qx,qval,phi,dx,dy,dist,two_theta 1453 1454 qval = CalcQval(xaxval,yaxval,xctr,yctr,sdd,lam,pixSize) 1455 1456 sdd *=100 //convert to cm 1457 dx = (xaxval  xctr)*pixSize //delta x in cm 1458 dy = (yaxval  yctr)*pixSize //delta y in cm 1459 phi = FindPhi(dx,dy) 1460 1461 //get scattering angle to project onto flat detector => Qr = qval*cos(theta) 1462 dist = sqrt(dx^2 + dy^2) 1463 two_theta = atan(dist/sdd) 1464 1465 qx = qval*cos(two_theta/2)*cos(phi) 1466 1467 return qx 1468 End 1469 1470 //calculates just the qvalue in the ydirection on the detector 1471 //input/output is the same as CalcQval() 1472 //ALL inputs are in detector coordinates 1473 //NOTE: detector locations passed in are pixel = 0.5cm real space on the Ordela detector 1474 //sdd is in meters 1475 //wavelength is in Angstroms 1476 // 1477 // repaired incorrect qx and qy calculation 3 dec 08 SRK (Lionel and C. Dewhurst) 1478 // now properly accounts for qz 1479 // 1480 Function CalcQY(xaxval,yaxval,xctr,yctr,sdd,lam,pixSize) 1481 Variable xaxval,yaxval,xctr,yctr,sdd,lam,pixSize 1482 1483 Variable dy,qval,dx,phi,qy,dist,two_theta 1484 1485 qval = CalcQval(xaxval,yaxval,xctr,yctr,sdd,lam,pixSize) 1486 1487 sdd *=100 //convert to cm 1488 dx = (xaxval  xctr)*pixSize //delta x in cm 1489 dy = (yaxval  yctr)*pixSize //delta y in cm 1490 phi = FindPhi(dx,dy) 1491 1492 //get scattering angle to project onto flat detector => Qr = qval*cos(theta) 1493 dist = sqrt(dx^2 + dy^2) 1494 two_theta = atan(dist/sdd) 1495 1496 qy = qval*cos(two_theta/2)*sin(phi) 1497 1498 return qy 1499 End 1500 1501 //calculates just the zcomponent of the qvector, not measured on the detector 1502 //input/output is the same as CalcQval() 1503 //ALL inputs are in detector coordinates 1504 //NOTE: detector locations passed in are pixel = 0.5cm real space on the Ordela detector 1505 //sdd is in meters 1506 //wavelength is in Angstroms 1507 // 1508 // not actually used, but here for completeness if anyone asks 1509 // 1510 Function CalcQZ(xaxval,yaxval,xctr,yctr,sdd,lam,pixSize) 1511 Variable xaxval,yaxval,xctr,yctr,sdd,lam,pixSize 1512 1513 Variable dy,qval,dx,phi,qz,dist,two_theta 1514 1515 qval = CalcQval(xaxval,yaxval,xctr,yctr,sdd,lam,pixSize) 1516 1517 sdd *=100 //convert to cm 1518 1519 //get scattering angle to project onto flat detector => Qr = qval*cos(theta) 1520 dx = (xaxval  xctr)*pixSize //delta x in cm 1521 dy = (yaxval  yctr)*pixSize //delta y in cm 1522 dist = sqrt(dx^2 + dy^2) 1523 two_theta = atan(dist/sdd) 1524 1525 qz = qval*sin(two_theta/2) 1526 1527 return qz 1528 End 1529 1530 //for commandline testing, replace the function declaration 1531 //Function FindQxQy(qq,phi) 1532 // Variable qq,phi 1533 // Variable qx,qy 1534 // 1535 // 1536 ThreadSafe Function FindQxQy(qq,phi,qx,qy) 1537 Variable qq,phi,&qx,&qy 1538 1539 qx = sqrt(qq^2/(1+tan(phi)*tan(phi))) 1540 qy = qx*tan(phi) 1541 1542 if(phi >= 0 && phi <= pi/2) 1543 qx = abs(qx) 1544 qy = abs(qy) 1545 endif 1546 1547 if(phi > pi/2 && phi <= pi) 1548 qx = abs(qx) 1549 qy = abs(qy) 1550 endif 1551 1552 if(phi > pi && phi <= pi*3/2) 1553 qx = abs(qx) 1554 qy = abs(qy) 1555 endif 1556 1557 if(phi > pi*3/2 && phi < 2*pi) 1558 qx = abs(qx) 1559 qy = abs(qy) 1560 endif 1561 1562 1563 // Print "recalculated qx,qy,q = ",qx,qy,sqrt(qx*qx+qy*qy) 1564 1565 return(0) 1566 end 1567
Note: See TracChangeset
for help on using the changeset viewer.