Changeset 634 for sans/XOP_Dev/SANSAnalysis/lib/libCylinder.c
 Timestamp:
 Mar 3, 2010 9:08:42 AM (13 years ago)
 File:

 1 edited
Legend:
 Unmodified
 Added
 Removed

sans/XOP_Dev/SANSAnalysis/lib/libCylinder.c
r632 r634 1919 1919 si1 = sin(sinarg1)/sinarg1; 1920 1920 } 1921 if ( besarg2 == 0.0){1921 if (sinarg2 == 0.0){ 1922 1922 si2 = 1.0; 1923 1923 }else{ … … 1968 1968 si1 = sin(sinarg1)/sinarg1; 1969 1969 } 1970 if ( besarg2 == 0.0){1970 if (sinarg2 == 0.0){ 1971 1971 si2 = 1.0; 1972 1972 }else{ … … 2058 2058 si1 = sin(sinarg1)/sinarg1; 2059 2059 } 2060 if ( besarg2 == 0.0){2060 if (sinarg2 == 0.0){ 2061 2061 si2 = 1.0; 2062 2062 }else{ … … 2456 2456 double summ,zi,yyy,answer; //running tally of integration 2457 2457 double summj,vaj,vbj,zij; //for the inner integration 2458 double SphCyl_tmp[7],arg1,arg2,inner ;2458 double SphCyl_tmp[7],arg1,arg2,inner,be; 2459 2459 2460 2460 … … 2506 2506 arg2 = x*rad*sin(zi); 2507 2507 yyy = inner; 2508 2509 if(arg2 == 0) { 2510 be = 0.5; 2511 } else { 2512 be = NR_BessJ1(arg2)/arg2; 2513 } 2508 2514 2509 2515 if(arg1 == 0.0) { //limiting value of sinc(0) is 1; sinc is not defined in math.h 2510 yyy += Pi*rad*rad*len*2.0* NR_BessJ1(arg2)/arg2;2516 yyy += Pi*rad*rad*len*2.0*be; 2511 2517 } else { 2512 yyy += Pi*rad*rad*len*sin(arg1)/arg1*2.0* NR_BessJ1(arg2)/arg2;2518 yyy += Pi*rad*rad*len*sin(arg1)/arg1*2.0*be; 2513 2519 } 2514 2520 yyy *= yyy; … … 2537 2543 double val,arg1,arg2; 2538 2544 double scale,bkg,sldc,slds; 2539 double len,rad,hDist,endRad ;2545 double len,rad,hDist,endRad,be; 2540 2546 scale = w[0]; 2541 2547 rad = w[1]; … … 2551 2557 arg2 = x*endRad*sin(theta)*sqrt(1.0tt*tt); 2552 2558 2553 val = cos(arg1)*(1.0tt*tt)*NR_BessJ1(arg2)/arg2; 2559 if(arg2 == 0) { 2560 be = 0.5; 2561 } else { 2562 be = NR_BessJ1(arg2)/arg2; 2563 } 2564 val = cos(arg1)*(1.0tt*tt)*be; 2554 2565 2555 2566 return(val); … … 2573 2584 double summ,zi,yyy,answer; //running tally of integration 2574 2585 double summj,vaj,vbj,zij; //for the inner integration 2575 double CLens_tmp[7],arg1,arg2,inner,hh ;2586 double CLens_tmp[7],arg1,arg2,inner,hh,be; 2576 2587 2577 2588 … … 2626 2637 yyy = inner; 2627 2638 2639 if(arg2 == 0) { 2640 be = 0.5; 2641 } else { 2642 be = NR_BessJ1(arg2)/arg2; 2643 } 2644 2628 2645 if(arg1 == 0.0) { //limiting value of sinc(0) is 1; sinc is not defined in math.h 2629 yyy += Pi*rad*rad*len*2.0* NR_BessJ1(arg2)/arg2;2646 yyy += Pi*rad*rad*len*2.0*be; 2630 2647 } else { 2631 yyy += Pi*rad*rad*len*sin(arg1)/arg1*2.0* NR_BessJ1(arg2)/arg2;2648 yyy += Pi*rad*rad*len*sin(arg1)/arg1*2.0*be; 2632 2649 } 2633 2650 yyy *= yyy; … … 2667 2684 double summ,zi,yyy,answer; //running tally of integration 2668 2685 double summj,vaj,vbj,zij; //for the inner integration 2669 double arg1,arg2,inner,hh ;2686 double arg1,arg2,inner,hh,be; 2670 2687 2671 2688 … … 2710 2727 yyy = inner; 2711 2728 2729 if(arg2 == 0) { 2730 be = 0.5; 2731 } else { 2732 be = NR_BessJ1(arg2)/arg2; 2733 } 2734 2712 2735 if(arg1 == 0.0) { //limiting value of sinc(0) is 1; sinc is not defined in math.h 2713 yyy += Pi*rad*rad*len*2.0* NR_BessJ1(arg2)/arg2;2736 yyy += Pi*rad*rad*len*2.0*be; 2714 2737 } else { 2715 yyy += Pi*rad*rad*len*sin(arg1)/arg1*2.0* NR_BessJ1(arg2)/arg2;2738 yyy += Pi*rad*rad*len*sin(arg1)/arg1*2.0*be; 2716 2739 } 2740 2741 2742 2717 2743 yyy *= yyy; 2718 2744 yyy *= sin(zi); // = A(q)^2*sin(theta) … … 2742 2768 double val,arg1,arg2; 2743 2769 double scale,bkg,sldc,slds; 2744 double len,rad,hDist,endRad ;2770 double len,rad,hDist,endRad,be; 2745 2771 scale = w[0]; 2746 2772 rad = w[1]; … … 2756 2782 arg2 = x*endRad*sin(theta)*sqrt(1.0tt*tt); 2757 2783 2758 val = cos(arg1)*(1.0tt*tt)*NR_BessJ1(arg2)/arg2; 2784 if(arg2 == 0) { 2785 be = 0.5; 2786 } else { 2787 be = NR_BessJ1(arg2)/arg2; 2788 } 2789 2790 val = cos(arg1)*(1.0tt*tt)*be; 2759 2791 2760 2792 return(val); … … 2778 2810 double summ,zi,yyy,answer; //running tally of integration 2779 2811 double summj,vaj,vbj,zij; //for the inner integration 2780 double Dumb_tmp[7],arg1,arg2,inner ;2812 double Dumb_tmp[7],arg1,arg2,inner,be; 2781 2813 2782 2814 … … 2831 2863 yyy = inner; 2832 2864 2865 if(arg2 == 0) { 2866 be = 0.5; 2867 } else { 2868 be = NR_BessJ1(arg2)/arg2; 2869 } 2870 2833 2871 if(arg1 == 0.0) { //limiting value of sinc(0) is 1; sinc is not defined in math.h 2834 yyy += Pi*rad*rad*len*2.0* NR_BessJ1(arg2)/arg2;2872 yyy += Pi*rad*rad*len*2.0*be; 2835 2873 } else { 2836 yyy += Pi*rad*rad*len*sin(arg1)/arg1*2.0* NR_BessJ1(arg2)/arg2;2874 yyy += Pi*rad*rad*len*sin(arg1)/arg1*2.0*be; 2837 2875 } 2838 2876 yyy *= yyy; … … 2872 2910 double summ,zi,yyy,answer; //running tally of integration 2873 2911 double summj,vaj,vbj,zij; //for the inner integration 2874 double arg1,arg2,inner ;2912 double arg1,arg2,inner,be; 2875 2913 2876 2914 … … 2915 2953 yyy = inner; 2916 2954 2955 if(arg2 == 0) { 2956 be = 0.5; 2957 } else { 2958 be = NR_BessJ1(arg2)/arg2; 2959 } 2960 2917 2961 if(arg1 == 0.0) { //limiting value of sinc(0) is 1; sinc is not defined in math.h 2918 yyy += Pi*rad*rad*len*2.0* NR_BessJ1(arg2)/arg2;2962 yyy += Pi*rad*rad*len*2.0*be; 2919 2963 } else { 2920 yyy += Pi*rad*rad*len*sin(arg1)/arg1*2.0* NR_BessJ1(arg2)/arg2;2964 yyy += Pi*rad*rad*len*sin(arg1)/arg1*2.0*be; 2921 2965 } 2922 2966 yyy *= yyy; … … 2948 2992 double val,arg1,arg2; 2949 2993 double scale,bkg,sldc,slds; 2950 double len,rad,hDist,endRad ;2994 double len,rad,hDist,endRad,be; 2951 2995 scale = w[0]; 2952 2996 rad = w[1]; … … 2962 3006 arg2 = x*endRad*sin(theta)*sqrt(1.0tt*tt); 2963 3007 2964 val = cos(arg1)*(1.0tt*tt)*NR_BessJ1(arg2)/arg2; 3008 if(arg2 == 0) { 3009 be = 0.5; 3010 } else { 3011 be = NR_BessJ1(arg2)/arg2; 3012 } 3013 val = cos(arg1)*(1.0tt*tt)*be; 2965 3014 2966 3015 return(val); … … 3059 3108 double sinarg1,sinarg2; 3060 3109 double t1,t2,t3; 3061 double retval ;3110 double retval,si1,si2,be1,be2; 3062 3111 3063 3112 double Pi = 4.0*atan(1.0); … … 3074 3123 sinarg2 = qq*(length+facthick)*cos(dum); 3075 3124 3076 t1 = 2.0*vol1*dr1*sin(sinarg1)/sinarg1*NR_BessJ1(besarg1)/besarg1; 3077 t2 = 2.0*vol2*dr2*sin(sinarg2)/sinarg2*NR_BessJ1(besarg2)/besarg2; 3078 t3 = 2.0*vol3*dr3*sin(sinarg2)/sinarg2*NR_BessJ1(besarg1)/besarg1; 3125 if(besarg1 == 0) { 3126 be1 = 0.5; 3127 } else { 3128 be1 = NR_BessJ1(besarg1)/besarg1; 3129 } 3130 if(besarg2 == 0) { 3131 be2 = 0.5; 3132 } else { 3133 be2 = NR_BessJ1(besarg2)/besarg2; 3134 } 3135 if(sinarg1 == 0) { 3136 si1 = 1.0; 3137 } else { 3138 si1 = sin(sinarg1)/sinarg1; 3139 } 3140 if(sinarg2 == 0) { 3141 si2 = 1.0; 3142 } else { 3143 si2 = sin(sinarg2)/sinarg2; 3144 } 3145 t1 = 2.0*vol1*dr1*si1*be1; 3146 t2 = 2.0*vol2*dr2*si2*be2; 3147 t3 = 2.0*vol3*dr3*si2*be1; 3079 3148 3080 3149 retval = ((t1+t2+t3)*(t1+t2+t3))*sin(dum);
Note: See TracChangeset
for help on using the changeset viewer.