- Timestamp:
- Nov 3, 2008 5:32:02 PM (14 years ago)
- Location:
- sans/Dev/trunk/NCNR_User_Procedures/Reduction/SANS
- Files:
-
- 2 edited
-
MultScatter_MonteCarlo_2D.ipf (modified) (7 diffs)
-
SASCALC.ipf (modified) (2 diffs)
Legend:
- Unmodified
- Added
- Removed
-
sans/Dev/trunk/NCNR_User_Procedures/Reduction/SANS/MultScatter_MonteCarlo_2D.ipf
r434 r436 12 12 // - Most importantly, this needs to be checked for correctness of the MC simulation 13 13 // X how can I get the "data" on absolute scale? This would be a great comparison vs. the ideal model calculation 14 // -why does my integrated tau not match up with John's analytical calculations? where are the assumptions?14 // X why does my integrated tau not match up with John's analytical calculations? where are the assumptions? 15 15 // - get rid of all small angle assumptions - to make sure that the calculation is correct at all angles 16 // - what is magical about Qu? Is this an assumpution? 17 18 // - at the larger angles, is the "flat" detector being properly accounted for - in terms of 16 17 // X at the larger angles, is the "flat" detector being properly accounted for - in terms of 19 18 // the solid angle and how many counts fall in that pixel. Am I implicitly defining a spherical detector 20 19 // so that what I see is already "corrected"? … … 27 26 // X if no MC desired, still use the selected model 28 27 // X better display of MC results on panel 29 // - settings for "count for X seconds" or "how long to 1E6 cts on detector" ( run short sim, then multiply)28 // - settings for "count for X seconds" or "how long to 1E6 cts on detector" (but 1E6 is typically too many counts...) 30 29 // - add quartz window scattering to the simulation somehow 31 30 // - do smeared models make any sense?? Yes, John agrees that they do, and may be used in a more realistic simulation … … 33 32 // or the volume fraction of solvent. 34 33 // 35 // -add to the results the fraction of coherently scattered neutrons that are singly scattered, different than34 // X add to the results the fraction of coherently scattered neutrons that are singly scattered, different than 36 35 // the overall fraction of singly scattered, and maybe more important to know. 37 36 // 38 // -change the fraction reaching the detector to exclude those that don't interact. These transmitted neutrons37 // X change the fraction reaching the detector to exclude those that don't interact. These transmitted neutrons 39 38 // aren't counted. Is the # that interact a better number? 40 39 // … … 160 159 results[6] = retWave[5]/retWave[1] //double coherent fraction 161 160 results[7] = retWave[6]/retWave[1] //multiple scatter fraction 162 results[8] = (retWave[0]-retWave[1])/retWave[0] //tra snmitted fraction161 results[8] = (retWave[0]-retWave[1])/retWave[0] //transmitted fraction 163 162 164 163 return(0) … … 219 218 results[6] = retWave[5]/retWave[1] //double coherent fraction 220 219 results[7] = retWave[6]/retWave[1] //multiple scatter fraction 221 results[8] = (retWave[0]-retWave[1])/retWave[0] //tra snmitted fraction220 results[8] = (retWave[0]-retWave[1])/retWave[0] //transmitted fraction 222 221 223 222 return(0) … … 339 338 Phi = 0.0 // Intialize azimuthal angle. 340 339 N1 += 1 // Increment total number neutrons counter. 341 DONE = 0 // True when neutron is absorbed or whenscattered out of the sample.340 DONE = 0 // True when neutron is scattered out of the sample. 342 341 INDEX = 0 // Set counter for number of scattering events. 343 342 zz = 0.0 // Set entering dimension of sample. … … 417 416 NN[INDEX] += 1 418 417 Endif 419 //IF (VZ > 1.0) // FIX INVALID ARGUMENT420 //VZ = 1.0 - 1.2e-7421 //ENDIF422 Theta_z = acos(Vz) // Angle WITH respect to z axis.423 testQ = 2*pi*sin(theta_z)/wavelength424 418 425 // pick a random phi angle, and see if it lands on the detector 426 // since the scattering is isotropic, I can safely pick a new, random value 427 // this would not be true if simulating anisotropic scattering. 428 testPhi = abs(enoise(1))*2*Pi 429 // is it on the detector? 430 FindPixel(testQ,testPhi,wavelength,sdd,pixSize,xCtr,yCtr,xPixel,yPixel) 419 if(index != 0) //the neutron interacted at least once, figure out where it ends up 420 421 Theta_z = acos(Vz) // Angle WITH respect to z axis. 422 testQ = 2*pi*sin(theta_z)/wavelength 423 424 // pick a random phi angle, and see if it lands on the detector 425 // since the scattering is isotropic, I can safely pick a new, random value 426 // this would not be true if simulating anisotropic scattering. 427 testPhi = abs(enoise(1))*2*Pi 428 // is it on the detector? 429 FindPixel(testQ,testPhi,wavelength,sdd,pixSize,xCtr,yCtr,xPixel,yPixel) 430 431 if(xPixel != -1 && yPixel != -1) 432 //if(index==1) // only the single scattering events 433 MC_linear_data[xPixel][yPixel] += 1 //this is the total scattering, including multiple scattering 434 //endif 435 isOn += 1 // neutron that lands on detector 436 endif 437 438 If(theta_z < theta_max) 439 //Choose index for scattering angle array. 440 //IND = NINT(THETA_z/DTH + 0.4999999) 441 ind = round(THETA_z/DTH + 0.4999999) //round is eqivalent to nint() 442 NT[ind] += 1 //Increment bin for angle. 443 //Increment angle array for single scattering events. 444 IF(INDEX == 1) 445 j1[ind] += 1 446 Endif 447 //Increment angle array for double scattering events. 448 IF (INDEX == 2) 449 j2[ind] += 1 450 Endif 451 EndIf 452 453 // increment all of the counters now since done==1 here and I'm sure to exit and get another neutron 454 NScatterEvents += index //total number of scattering events 455 if(index == 1 && incoherentEvent == 1) 456 NSingleIncoherent += 1 457 endif 458 if(index == 1 && coherentEvent == 1) 459 NSingleCoherent += 1 460 endif 461 if(index == 2 && coherentEvent == 1 && incoherentEvent == 0) 462 NDoubleCoherent += 1 463 endif 464 if(index > 1) 465 NMultipleScatter += 1 466 endif 467 //Print "n1,index (x,y) = ",n1,index, xpixel,ypixel 468 else // if neutron escaped without interacting 469 // then it must be a transmitted neutron 470 // don't need to calculate, just increment the proper counters 471 MC_linear_data[xCtr][yCtr] += 1 472 isOn += 1 473 nt[0] += 1 474 475 // (Don't do this - it's a waste of time) 476 // re-initialize and go back and count this neutron again 477 // go back to xy position 478 // xx -= ll*vx 479 // yy -= ll*vy 480 // zz -= ll*vz 481 // RR = sqrt(xx*xx+yy*yy) //radial position of scattering event. 482 // 483 // Vx = 0.0 // Initialize incident direction vector. 484 // Vy = 0.0 485 // Vz = 1.0 486 // zz = 0 487 // theta = 0 488 // phi = 0 489 endif 431 490 432 if(xPixel != -1 && yPixel != -1)433 //if(index==1) // only the single scattering events434 MC_linear_data[xPixel][yPixel] += 1 //this is the total scattering, including multiple scattering435 //endif436 isOn += 1 // scattered neutron that lands on detector437 endif438 439 If(theta_z < theta_max)440 //Choose index for scattering angle array.441 //IND = NINT(THETA_z/DTH + 0.4999999)442 ind = round(THETA_z/DTH + 0.4999999) //round is eqivalent to nint()443 NT[ind] += 1 //Increment bin for angle.444 //Increment angle array for single scattering events.445 IF(INDEX == 1)446 j1[ind] += 1447 Endif448 //Increment angle array for double scattering events.449 IF (INDEX == 2)450 j2[ind] += 1451 Endif452 EndIf453 454 // increment all of the counters now since done==1 here and I'm sure to exit and get another neutron455 NScatterEvents += index //total number of scattering events456 if(index == 1 && incoherentEvent == 1)457 NSingleIncoherent += 1458 endif459 if(index == 1 && coherentEvent == 1)460 NSingleCoherent += 1461 endif462 if(index == 2 && coherentEvent == 1 && incoherentEvent == 0)463 NDoubleCoherent += 1464 endif465 if(index > 1)466 NMultipleScatter += 1467 endif468 491 ENDIF 469 492 while (!done) -
sans/Dev/trunk/NCNR_User_Procedures/Reduction/SANS/SASCALC.ipf
r434 r436 711 711 inputWave[10] = sig_sas 712 712 713 linear_data = 0 //initialize 714 713 715 Variable t0 = stopMStimer(-2) 714 716 … … 722 724 Variable trans 723 725 trans = results[8] //(n1-n2)/n1 726 if(trans == 0) 727 trans = 1 728 endif 724 729 725 730 // convert to absolute scale
Note: See TracChangeset
for help on using the changeset viewer.
