Estimate transmission wrapper#365
Conversation
| transmission_limits = (get_beamline_param( | ||
| beamline_config, ("gda.mx.udc.minTransmission",), 0.0), | ||
| min( | ||
| get_beamline_param( | ||
| beamline_config, ("gda.mx.udc.maxTransmission",), 1.0 | ||
| ), | ||
| recommended_max_transmission | ||
| ), | ||
| ) |
There was a problem hiding this comment.
Discussed on slack that rather than using the recommended_max transmission to set the transmission limits, the recommended transmission should be used in place of recipe_step.transmission (i.e. instead of the transmission value from the Agamemnon recipe). This means that the recommended transmission will be scaled appropriately for wavelength and resolution.
There was a problem hiding this comment.
I have tested with Phasing which has a recommended transmission of .25, and the recommended transmission from our wrapper is .44 and after scaling it returns to use .36 transmission in the dc. You can see it at https://ispyb-test.diamond.ac.uk/dc/visit/mx23694-156/id/21149301 for Phasing wedge 1
| for command, script, output_file in commands: | ||
| result = subprocess.run(script, cwd=working_directory, check=True) | ||
|
|
||
| while not os.path.exists(output_file): |
There was a problem hiding this comment.
Is this step necessary? Have you found that the dials.import subprocess returns without the output files existing?
|
|
||
| reflection_file = results_directory / "strong.refl" | ||
| reflections = flex.reflection_table.from_file(reflection_file) | ||
| num_counts, _ = self.build_hist_from_reflections(reflections) |
There was a problem hiding this comment.
As you've gone to the trouble of constructing a histogram, it would be good to plot this and have it as an output file. Additionally, you could output the histogram similarly to how xia2.overload does, so that the data is available as a result file. If this takes a some time to process, you could always run this after the wrapper sends the scaled_transmission off to the next recipe step, so that the pipeline is not slowed down.
There was a problem hiding this comment.
The histogram can just be saved using json.dump(...) and does not take much time. Would there be another why using JSONLines is preferable?
There was a problem hiding this comment.
No reason to use JSONLines, you can just create it in the wrapper as you have done
| return False | ||
|
|
||
| results_directory.mkdir(parents=True, exist_ok=True) | ||
| output_files = ["dials.find_spots.log", "strong.refl", "imported.expt"] |
There was a problem hiding this comment.
The strong.refl and imported.expt are probably not too useful to have as results files so I would omit them from here. I think it would be useful to construct a log file for this wrapper, that gives a simple overview of the running of the dials commands as well as the analysis of the spots afterwards.
| curr_program = find_process_program.filter( | ||
| AutoProcProgram.autoProcProgramId == parameters.program_id | ||
| ).scalar() | ||
| # xia2 dials occassionaly gives optimistic estimate for resolution | ||
| if curr_program == "xia2 dials": | ||
| self.log.info( | ||
| f"Skipping strategy trigger for dcid={parameters.dcid} from program: xia2 dials." | ||
| ) | ||
| return {"success": True} |
There was a problem hiding this comment.
This query/check isn't needed. We can avoid triggering for xia2-dials by simply not including the call to the trigger service in the xia2-dials recipe(s).
| def draw_plot(self, counts, pixels): | ||
| """Create ASCII art histogram""" | ||
|
|
||
| self.log.info("Plotting pixel intensities...") | ||
|
|
||
| width, height = 60, 20 | ||
| title = "'Pixel intensity distribution'" | ||
| xlabel = "'Num counts'" | ||
| ylabel = "'Counts'" | ||
|
|
||
| command = ["gnuplot"] | ||
| plot_commands = [ | ||
| "set term dumb %d %d" % (width, height - 2), | ||
| "set logscale y", | ||
| "set logscale x", | ||
| "set ytics out", | ||
| "set title %s" % title, | ||
| "set xlabel %s" % xlabel, | ||
| "set ylabel %s offset character %d,0" % (ylabel, len(ylabel) // 2), | ||
| ] | ||
|
|
||
| data_string = "\n".join(f"{x} {y}" for x, y in zip(counts[1:], pixels[1:])) | ||
|
|
||
| plot_commands.append("plot '-' using 1:2 title '' with lines") | ||
| plot_commands.append(data_string) | ||
|
|
||
| try: | ||
| result = subprocess.run( | ||
| command, input="\n".join(plot_commands), text=True, capture_output=True | ||
| ) | ||
| except (OSError, subprocess.TimeoutExpired) as e: | ||
| self.log.info("Error plotting counts vs pixels") | ||
| self.log.info(e) | ||
| return | ||
| else: | ||
| self.log.info(result.stdout) |
There was a problem hiding this comment.
I'd like to avoid having gnuplot as a dependency. I'd be happy just to have a matplotlib png of the histogram as an output.
| return len(num_pixels) | ||
|
|
||
| def save_hist_to_json(self, hist, max_trusted_value, results_dir): | ||
| results_path = results_dir / "overload.json" |
There was a problem hiding this comment.
We could give this a more informative name and don't need to follow the xia2_overload filename. Something like pixel_count_histogram.json
3 participants
Wrapper to estimate the max recommended transmission from dials.find_spots.