Overall the two reviewers provided positive comments on the quality and novelty of the manuscript, and the authors provided diligent responses to each of the authors' comments. I think the authors have adequately addressed most of the reviewers' concerns, except one concern of Reviewer 2, that is, they authors only tested the new method in one watershed. I share the same concern, therefore hope the authors will add one more watershed in the manuscript before submitting back the revised manuscript for further review.

The authors have revised the manuscript by adding more test watersheds as suggested by reviewer 2. As the reviewer has major concerns over this issue, I will send it for the reviewers to conduct the second round review.

The two reviewers gave contrasting comments about the quality of the manuscript. Reviewer 1, while appreciate the responses of the authors, still have serious concerns about the methods used and the procedures used to compare different methods. By balancing the two reviewers comments, I think there is still new value in the manuscript, but also agree with reviewer 1 that additional substantial revisions are necessary in order to bring the manuscript to meet the publication bar. My specific comments are as follows:

[1] Regarding reviewer 1’s first comment “Response to comments 29: you mention that the hydroclimatological model requires spatially varied parameters, but you only used one set "to verify the model with less calibration effort" in Athabasca. I wonder what the results would be if you calibrated the hydroloclimatological model as it should have been done. Would the results be similar to the North Fork Tolt River?”:

Please confirm if you used the same calibration procedures to calibrate the models in the two watersheds. If they are not the same, the authors must provide strong justification why they use different procedures.

[2] Regarding reviewer 1’s second comment “Response to comments 32: I do not agree with your explanation. You write that "lower water depth...means stream temperature changes quicker". In SWAT, river width can be orders of magnitude off the real value because the simplified GIS algorithm does not account for local stream characteristics. You can easily check this by overlaying your stream network on Satellite Images. In most cases this is not a significant problem (though it can be for routing and transmission losses). But in your case, if your width is overestimated by a factor of let's say 5, your depth would be underestimated by a similar value. Keeping in mind that you are proposing a new method of stream temperature calculations in SWAT, I think you should know
(a) if your river width is somewhere near reality and if not, correct it in the .rte (which is not complicated: you can cluster your subbasins according to river order and get representative widths for each order and change the .rte accordingly) and
(b) how sensitivity water depth is in your equation”:

Please check if the bankful width used the model match reality values. If credible data is not available, the authors should discuss potential uncertainties associated with the algorithms used in SWAT. I understand that correcting the errors associated with the built-in bankful width estimation algorithm within SWAT is beyond the scope of the study, but certain sensitivity test would be helpful to understand the value river width information for accurate water temperature simulation.

[3] Regarding reviewer 1’s third comment “Response to Comment 38: Sorry, your answer and p.9 l.16-19 is still confusing for me. "blue dots are the observed period daily average (not continuous) ... measured periodically (monthly or weekly)" So how long is an "observed period daily average" continuously measured? I don't think someone is waiting at each gauge for 24hours and takes a sample every hour, but I suspect what you got is 'only' an instantaneous observation at the particular time of the day when the observation took place at max once a month, or? In that case, I think your explanation is misleading. Please rewrite so that the reader knows the true quality of your data. Considering this, I do not agree with your abstract p.1 l.21-22: What you used is not "high frequency observed stream temperature data"”:

Please better clarify the description of the data used, and change the abstract to reflect the actual quality of the data you used.

[4] Regarding reviewer 2’s fourth comment “p.14 l.24 - p.15 l.9: Your model performs worse in the second catchment. This raises the question "Why is this the case"? You don't discuss possible reasons anywhere.
And equally important: why did you choose this particular catchment out of the seven other catchments? Do you expect your model to perform better in the other catchments?
Besides that, in the abstract and conclusion you only mention the RE and hide the NSE, but in the Athabasca where the NSE is higher, you mention the NSE. I think this is not good practice.”:

Please clarify the performance of different methods. I also found it is difficult to understand the comparative analysis of the performances of different methods. One suggestion is to use “default” to represent the SWAT built-in method, use “Ficklin et al.’s method” to indicate the work by Ficklin et al., and use “this study” to represent the new method proposed in the manuscript.

If NSE or RE is worse for the new method, the authors need to identify strong reasons to justify why the new method is better than the existing ones. If the new method performed less in the second catchment, and the authors would like justify the new method’s advantage to the existing ones, more watershed should be used and the new method should perform better for most cases.

[5] Regarding reviewer 1’s figth comment “p.17 Figure 6: It looks like the original SWAT is exactly below the newly developed model, right? Is this a labelling error? Considering that you built the model on top of the Hydroclimatological model (which is significantly different from the two), this is very confusing. Based on this result, I do not understand your interpretation on p.18, l.8-9 "the equilibrium temp provides a potential tool for more-accurate water quality concentration simulation" Based on these results, I think it is impossible to know which model version is "more-accurate".
Also, the selection of lines and markers in unfortunate. I suggest to remove the markers and use different thicknesses for the two lines that closely ovelap.”:

I think the purpose of figure 6 is to show that organic N simulation is sensitive to temperature module. Anyway, please add more information in the caption to highlight the key message from this figure. Are they used to show the new method is better, or the authors intend to show the sensitivity of organic N simulation? Also, please improve the legibility, and clearly label the names of the different method, and which one is the method proposed in the manuscript.