@@ -712,7 +712,7 @@ void Monitor_nd_noaccInit(Monitornd_noaccDefines_type *DEFS,
712712 { /* Dim : (Vars->Coord_Number+1)*Vars->Buffer_Block matrix (for p, dp) */
713713 Vars -> Mon2D_Buffer = (double * )malloc ((Vars -> Coord_Number + 1 )* Vars -> Buffer_Block * sizeof (double ));
714714 if (Vars -> Mon2D_Buffer == NULL )
715- { printf ("Monitor_nD: %s cannot allocate Vars->Mon2D_Buffer (%li ). No list and auto limits.\n" , Vars -> compcurname , Vars -> Buffer_Block * (Vars -> Coord_Number + 1 )* sizeof (double )); Vars -> Flag_List = 0 ; Vars -> Flag_Auto_Limits = 0 ; }
715+ { printf ("Monitor_nD: %s cannot allocate Vars->Mon2D_Buffer (%zi ). No list and auto limits.\n" , Vars -> compcurname , Vars -> Buffer_Block * (Vars -> Coord_Number + 1 )* sizeof (double )); Vars -> Flag_List = 0 ; Vars -> Flag_Auto_Limits = 0 ; }
716716 else
717717 {
718718 for (i = 0 ; i < (Vars -> Coord_Number + 1 )* Vars -> Buffer_Block ; Vars -> Mon2D_Buffer [i ++ ] = (double )0 );
@@ -727,14 +727,14 @@ void Monitor_nd_noaccInit(Monitornd_noaccDefines_type *DEFS,
727727 Vars -> Mon2D_p = (double * * )malloc ((Vars -> Coord_Number )* sizeof (double * ));
728728 Vars -> Mon2D_p2 = (double * * )malloc ((Vars -> Coord_Number )* sizeof (double * ));
729729 if ((Vars -> Mon2D_N == NULL ) || (Vars -> Mon2D_p == NULL ) || (Vars -> Mon2D_p2 == NULL ))
730- { fprintf (stderr ,"Monitor_nD: %s n1D cannot allocate Vars->Mon2D_N/p/p2 (%li ). Fatal.\n" , Vars -> compcurname , (Vars -> Coord_Number )* sizeof (double * )); exit (-1 ); }
730+ { fprintf (stderr ,"Monitor_nD: %s n1D cannot allocate Vars->Mon2D_N/p/p2 (%zi ). Fatal.\n" , Vars -> compcurname , (Vars -> Coord_Number )* sizeof (double * )); exit (-1 ); }
731731 for (i = 1 ; i <= Vars -> Coord_Number ; i ++ )
732732 {
733733 Vars -> Mon2D_N [i - 1 ] = (double * )malloc (Vars -> Coord_Bin [i ]* sizeof (double ));
734734 Vars -> Mon2D_p [i - 1 ] = (double * )malloc (Vars -> Coord_Bin [i ]* sizeof (double ));
735735 Vars -> Mon2D_p2 [i - 1 ] = (double * )malloc (Vars -> Coord_Bin [i ]* sizeof (double ));
736736 if ((Vars -> Mon2D_N == NULL ) || (Vars -> Mon2D_p == NULL ) || (Vars -> Mon2D_p2 == NULL ))
737- { fprintf (stderr ,"Monitor_nD: %s n1D cannot allocate %s Vars->Mon2D_N/p/p2[%li] (%li ). Fatal.\n" , Vars -> compcurname , Vars -> Coord_Var [i ], i , (Vars -> Coord_Bin [i ])* sizeof (double * )); exit (-1 ); }
737+ { fprintf (stderr ,"Monitor_nD: %s n1D cannot allocate %s Vars->Mon2D_N/p/p2[%li] (%zi ). Fatal.\n" , Vars -> compcurname , Vars -> Coord_Var [i ], i , (Vars -> Coord_Bin [i ])* sizeof (double * )); exit (-1 ); }
738738 else
739739 {
740740 for (j = 0 ; j < Vars -> Coord_Bin [i ]; j ++ )
@@ -749,14 +749,14 @@ void Monitor_nd_noaccInit(Monitornd_noaccDefines_type *DEFS,
749749 Vars -> Mon2D_p = (double * * )malloc ((Vars -> Coord_Bin [1 ])* sizeof (double * ));
750750 Vars -> Mon2D_p2 = (double * * )malloc ((Vars -> Coord_Bin [1 ])* sizeof (double * ));
751751 if ((Vars -> Mon2D_N == NULL ) || (Vars -> Mon2D_p == NULL ) || (Vars -> Mon2D_p2 == NULL ))
752- { fprintf (stderr ,"Monitor_nD: %s 2D cannot allocate %s Vars->Mon2D_N/p/p2 (%li ). Fatal.\n" , Vars -> compcurname , Vars -> Coord_Var [1 ], (Vars -> Coord_Bin [1 ])* sizeof (double * )); exit (-1 ); }
752+ { fprintf (stderr ,"Monitor_nD: %s 2D cannot allocate %s Vars->Mon2D_N/p/p2 (%zi ). Fatal.\n" , Vars -> compcurname , Vars -> Coord_Var [1 ], (Vars -> Coord_Bin [1 ])* sizeof (double * )); exit (-1 ); }
753753 for (i = 0 ; i < Vars -> Coord_Bin [1 ]; i ++ )
754754 {
755755 Vars -> Mon2D_N [i ] = (double * )malloc (Vars -> Coord_Bin [2 ]* sizeof (double ));
756756 Vars -> Mon2D_p [i ] = (double * )malloc (Vars -> Coord_Bin [2 ]* sizeof (double ));
757757 Vars -> Mon2D_p2 [i ] = (double * )malloc (Vars -> Coord_Bin [2 ]* sizeof (double ));
758758 if ((Vars -> Mon2D_N == NULL ) || (Vars -> Mon2D_p == NULL ) || (Vars -> Mon2D_p2 == NULL ))
759- { fprintf (stderr ,"Monitor_nD: %s 2D cannot allocate %s Vars->Mon2D_N/p/p2[%li] (%li ). Fatal.\n" , Vars -> compcurname , Vars -> Coord_Var [1 ], i , (Vars -> Coord_Bin [2 ])* sizeof (double * )); exit (-1 ); }
759+ { fprintf (stderr ,"Monitor_nD: %s 2D cannot allocate %s Vars->Mon2D_N/p/p2[%li] (%zi ). Fatal.\n" , Vars -> compcurname , Vars -> Coord_Var [1 ], i , (Vars -> Coord_Bin [2 ])* sizeof (double * )); exit (-1 ); }
760760 else
761761 {
762762 for (j = 0 ; j < Vars -> Coord_Bin [2 ]; j ++ )
@@ -810,6 +810,82 @@ void Monitor_nd_noaccInit(Monitornd_noaccDefines_type *DEFS,
810810 if (mpi_node_rank == mpi_node_root ) {
811811 #endif
812812 if (nxhandle ) {
813+
814+ /* This section of code writes detector shape information to
815+ entryN/instrument/components/'name'/geometry in the NeXus file */
816+
817+ char nexuscomp [CHAR_BUF_LENGTH ];
818+ char pref [5 ];
819+ if (Vars -> compcurindex - 1 < 10 ) {
820+ sprintf (pref ,"000" );
821+ } else if (Vars -> compcurindex - 1 < 100 ) {
822+ sprintf (pref ,"00" );
823+ } else if (Vars -> compcurindex - 1 < 1000 ) {
824+ sprintf (pref ,"0" );
825+ } else if (Vars -> compcurindex - 1 < 10000 ) {
826+ sprintf (pref ,"" );
827+ } else {
828+ fprintf (stderr ,"Error, no support for > 10000 comps at the moment!\n" );
829+ exit (-1 );
830+ }
831+ sprintf (nexuscomp ,"%s%d_%s" ,pref ,Vars -> compcurindex - 1 ,Vars -> compcurname );
832+
833+ if (NXopengroup (nxhandle , "instrument" , "NXinstrument" ) == NX_OK ) {
834+ if (NXopengroup (nxhandle , "components" , "NXdata" ) == NX_OK ) {
835+ if (NXopengroup (nxhandle , nexuscomp , "NXdata" ) == NX_OK ) {
836+ if (NXmakegroup (nxhandle , "Geometry" , "NXdata" ) == NX_OK ) {
837+ if (NXopengroup (nxhandle , "Geometry" , "NXdata" ) == NX_OK ) {
838+ char tmp [CHAR_BUF_LENGTH ];
839+ sprintf (tmp ,"%g" ,Vars -> Sphere_Radius );
840+ nxprintattr (nxhandle , "radius" , tmp );
841+
842+ sprintf (tmp ,"%g" ,Vars -> Cylinder_Height );
843+ nxprintattr (nxhandle , "height" , tmp );
844+
845+ sprintf (tmp ,"%g" ,Vars -> mxmin );
846+ nxprintattr (nxhandle , "xmin" , tmp );
847+ sprintf (tmp ,"%g" ,Vars -> mxmax );
848+ nxprintattr (nxhandle , "xmax" , tmp );
849+
850+ sprintf (tmp ,"%g" ,Vars -> mymin );
851+ nxprintattr (nxhandle , "ymin" , tmp );
852+ sprintf (tmp ,"%g" ,Vars -> mymax );
853+ nxprintattr (nxhandle , "ymax" , tmp );
854+
855+ sprintf (tmp ,"%g" ,Vars -> mzmin );
856+ nxprintattr (nxhandle , "zmin" , tmp );
857+ sprintf (tmp ,"%g" ,Vars -> mzmax );
858+ nxprintattr (nxhandle , "zmax" , tmp );
859+
860+ sprintf (tmp ,"%g" ,Vars -> mzmin );
861+ nxprintattr (nxhandle , "zmin" , tmp );
862+ sprintf (tmp ,"%g" ,Vars -> mzmax );
863+ nxprintattr (nxhandle , "zmax" , tmp );
864+
865+ sprintf (tmp ,"%i" ,Vars -> Flag_Shape );
866+ nxprintattr (nxhandle , "Shape identifier" , tmp );
867+ sprintf (tmp ,"%s" ,Vars -> Monitor_Label );
868+ nxprintattr (nxhandle , "Shape string" , tmp );
869+ sprintf (tmp ,"%s" ,Vars -> option );
870+ nxprintattr (nxhandle , "Option string" , tmp );
871+
872+ NXclosegroup (nxhandle ); // Geometry
873+ } else {
874+ printf ("Failed to open component NeXus component Geometry group\n" );
875+ }
876+ } else {
877+ printf ("Failed to create component NeXus component Geometry group\n" );
878+ }
879+ NXclosegroup (nxhandle ); // component
880+ }
881+ NXclosegroup (nxhandle ); // components
882+ } else {
883+ printf ("Failed to open NeXus component hierarchy\n" );
884+ }
885+ NXclosegroup (nxhandle ); // instrument
886+ }
887+
888+ /* Below code communicates geometry-oriented "BINS" for the detector. */
813889 char metadata [CHAR_BUF_LENGTH ];
814890 char metadatatmp [CHAR_BUF_LENGTH ];
815891 // Vars for 1D, >3D, OFF
@@ -820,9 +896,17 @@ void Monitor_nd_noaccInit(Monitornd_noaccDefines_type *DEFS,
820896 char binvar [CHAR_BUF_LENGTH ];
821897 sprintf (binlabel ,"none" );
822898 sprintf (binvar ,"none" );
823- long pix = Vars -> Coord_Min [Vars -> Coord_Number - 1 ]; // Second to last col is min. pixel id
824899
900+ // Find index of pixel column
901+ int id_index ;
902+ for (id_index = 0 ;id_index < 30 ;id_index ++ ) {
903+ if (strcmp (Vars -> Coord_Var [id_index ], "id" ) == 0 ) break ;
904+ }
905+ if (id_index == 30 ) id_index = Vars -> Coord_Number - 1 ; // Revert to earlier behavior is id not found
906+ long pix = Vars -> Coord_Min [id_index ];
907+
825908 MCDETECTOR detector ;
909+
826910 /* Init - perhaps better with an init-function in mccode-r? */
827911 detector .m = 0 ;
828912 detector .xmin = 0 ;
@@ -857,7 +941,7 @@ void Monitor_nd_noaccInit(Monitornd_noaccDefines_type *DEFS,
857941
858942 if (!Vars -> Flag_OFF ) {
859943
860- sprintf (metadata ,"id=%ld + %ld pixels: " ,(long )Vars -> Coord_Min [Vars -> Coord_Number - 1 ],(long )Vars -> Coord_BinProd [Vars -> Coord_Number ]);
944+ sprintf (metadata ,"id=%ld + %ld pixels: " ,(long )Vars -> Coord_Min [id_index ],(long )Vars -> Coord_BinProd [Vars -> Coord_Number ]);
861945 for (i = 1 ; i < N_spatial_dims ; i ++ ) {
862946 sprintf (metadatatmp ,"%s %s (%ld bins) x " ,metadata ,Vars -> Coord_Label [i ],Vars -> Coord_Bin [i ]);
863947 sprintf (metadata ,"%s" ,metadatatmp );
@@ -999,9 +1083,6 @@ int Monitor_nd_noaccTrace(Monitornd_noaccDefines_type *DEFS, Monitornd_noaccVari
9991083 char While_End = 0 ;
10001084 long While_Buffer = 0 ;
10011085 char Set_Vars_Coord_Type = DEFS -> COORD_NONE ;
1002-
1003- long long ParticleCount ;
1004- ParticleCount = Vars -> Neutron_Counter ++ ;
10051086
10061087 /* the logic below depends mainly on:
10071088 Flag_List: 1=store 1 buffer, 2=list all, 3=re-use buffer
@@ -1051,7 +1132,10 @@ int Monitor_nd_noaccTrace(Monitornd_noaccDefines_type *DEFS, Monitornd_noaccVari
10511132 Vars -> Mon2D_Buffer = (double * )realloc (Vars -> Mon2D_Buffer , (Vars -> Coord_Number + 1 )* (ParticleCount + Vars -> Buffer_Block )* sizeof (double ));
10521133 if (Vars -> Mon2D_Buffer == NULL )
10531134 { printf ("Monitor_nD: %s cannot reallocate Vars->Mon2D_Buffer[%li] (%li). Skipping.\n" , Vars -> compcurname , i , (long int )(ParticleCount + Vars -> Buffer_Block )* sizeof (double )); Vars -> Flag_List = 1 ; }
1054- else { Vars -> Buffer_Counter = 0 ; Vars -> Buffer_Size = ParticleCount + Vars -> Buffer_Block ; }
1135+ else {
1136+ Vars -> Buffer_Block = 2 * Vars -> Buffer_Block ;
1137+ Vars -> Buffer_Size = Vars -> Buffer_Block ;
1138+ }
10551139 }
10561140 } /* end if Buffer realloc */
10571141#endif
@@ -1102,7 +1186,8 @@ int Monitor_nd_noaccTrace(Monitornd_noaccDefines_type *DEFS, Monitornd_noaccVari
11021186 }
11031187 if (!flag_outside ) {
11041188 Vars -> Mon2D_Buffer [i + While_Buffer * (Vars -> Coord_Number + 1 )] = Coord [i ];
1105- pixID = Coord [i ];
1189+ printf ("Setting pixel ID to %g\n" ,Coord [i ]);
1190+ * pixID = Coord [i ];
11061191 }
11071192 } /* end if PixelID */
11081193 }
@@ -1294,7 +1379,9 @@ int Monitor_nd_noaccTrace(Monitornd_noaccDefines_type *DEFS, Monitornd_noaccVari
12941379 XY += Coord_Index [j ]* Vars -> Coord_BinProd [j - 1 ];
12951380 }
12961381 if (Vars -> Flag_mantid && Vars -> Flag_OFF && Vars -> OFF_polyidx >=0 ) XY = Vars -> OFF_polyidx ;
1297- if (!flag_outside ) XY += Vars -> Coord_Min [i ];
1382+ if (!flag_outside ) {
1383+ XY += Vars -> Coord_Min [i ];
1384+ }
12981385 }
12991386
13001387 /* handle 'abs' and 'log' keywords */
@@ -1386,9 +1473,8 @@ int Monitor_nd_noaccTrace(Monitornd_noaccDefines_type *DEFS, Monitornd_noaccVari
13861473 {
13871474 for (i = 0 ; i <= Vars -> Coord_Number ; i ++ )
13881475 {
1389- Vars -> Mon2D_Buffer [i + ParticleCount * (Vars -> Coord_Number + 1 )] = Coord [i ];
1476+ Vars -> Mon2D_Buffer [i + Vars -> Buffer_Counter * (Vars -> Coord_Number + 1 )] = Coord [i ];
13901477 }
1391- #pragma acc atomic
13921478 Vars -> Buffer_Counter = Vars -> Buffer_Counter + 1 ;
13931479 if (Vars -> Flag_Verbose && (Vars -> Buffer_Counter >= Vars -> Buffer_Block ) && (Vars -> Flag_List == 1 ))
13941480 printf ("Monitor_nD: %s %li neutrons stored in List.\n" , Vars -> compcurname , Vars -> Buffer_Counter );
@@ -1410,6 +1496,8 @@ int Monitor_nd_noaccTrace(Monitornd_noaccDefines_type *DEFS, Monitornd_noaccVari
14101496 }else {
14111497 return -1 ;
14121498 }
1499+ } else {
1500+
14131501 }
14141502 return 1 ;
14151503} /* end Monitor_nd_noaccTrace */
@@ -1716,7 +1804,7 @@ MCDETECTOR Monitor_nd_noaccSave(Monitornd_noaccDefines_type *DEFS, Monitornd_noa
17161804 p2m = (double * )malloc (Vars -> Coord_Bin [1 ]* Vars -> Coord_Bin [2 ]* sizeof (double ));
17171805 if (p2m == NULL )
17181806 {
1719- if (Vars -> Flag_Verbose ) printf ("Monitor_nD: %s cannot allocate memory for 2D array (%li ). Skipping.\n" , Vars -> compcurname , 3 * Vars -> Coord_Bin [1 ]* Vars -> Coord_Bin [2 ]* sizeof (double ));
1807+ if (Vars -> Flag_Verbose ) printf ("Monitor_nD: %s cannot allocate memory for 2D array (%zi ). Skipping.\n" , Vars -> compcurname , 3 * Vars -> Coord_Bin [1 ]* Vars -> Coord_Bin [2 ]* sizeof (double ));
17201808 /* comment out 'free memory' lines to avoid loosing arrays if
17211809 'detector' structure is used by other instrument parts
17221810 if (p0m != NULL) free(p0m);
0 commit comments